Implementation of Aerosol-Cloud Interaction within WRF-CHIMERE Online Coupled Model: Evaluation and Investigation of the Indirect Radiative Effect from Anthropogenic Emission Reduction on the Benelux Union

of Aerosol-Cloud Interaction within WRF-CHIMERE Online Investigation the Indirect Effect Anthropogenic Abstract: The indirect effects of aerosol are particularly important over regions where meteorological conditions and aerosol content are favourable to cloud formation. This was observed during the Intensive Cloud Aerosol Measurement Campaign (IMPACT) (European Integrated project on Aerosol Cloud Climate and Air quality Interaction (EUCAARI) project) in the Benelux Union during May 2008. To better understand this cloud formation variability, the indirect effects of aerosol have been included within the WRF-CHIMERE online model. By comparing model results to the aircraft measurements of IMPACT, to surface measurements from EMEP and AIRBASE and to MODIS satellite measurements, we showed that the model is able to simulate the variability and order of magnitude of the observed number of condensation nuclei (CN), even if some differences are identiﬁed for speciﬁc aerosol size and location. To quantify the impact of the local anthropogenic emissions on cloud formation, a sensitivity study is performed by halving the surface emissions ﬂuxes. It is shown that the indirect radiative effect (IRE) at the surface is positive for both shortwave and longwave with a net warming of +0.99 W/m 2 . In addition, important instantaneous changes are modelled at local scale with up to ± 6 ◦ C for temperatures and ± 50 mm/day for precipitation. R.B. and L.M.; Performed model simulations: P.T.; Software P.T., R.B., L.M., S.M., G.S. and D.K.; Methodology and model evaluation, P.T.; Formal analysis, P.T.; writing—original draft preparation, P.T.; writing—review and editing, P.T., L.M., A.D., S.M., S.T. and R.B.; supervision, L.M.; project administration, L.M.; funding acquisition, L.M.

determined by using the parameterization of DeMott et al. [47] that is a modified version of the scheme of De Mott et al. [11]. The freezing of supercooled water by contact nucleation follows Bigg [48] as done by Thompson et al. [49] with the IN concentration that alter the "effective temperature" to freeze the water drops. Homogeneous freezing of deliquesced aerosols is parameterized following Koop et al. [50]. Deliquesced aerosol number concentration is based on the CHIMERE prediction of a mixture of hygroscopic particles (sulphate, nitrate, ammonium, sea salt and secondary organic aerosols) with the diameter larger than 0.1 μm.
In summary, size number distribution, aerosol bulk hygroscopitcity, IN and deliquesced aerosol predicted by CHIMERE are sent to WRF through the coupler OASIS. These fields are used within WRF to calculate the droplet and ice nucleation rates. Cloud droplet and ice number concentration affect both radiative and microphysical processes. Radiative transfer is forced by cloud optical depth calculated starting from predicted cloud droplet and ice effective radius. Moreover, cloud microphysical processes rates are affected by cloud droplet and ice number concentration. Finally, the resulting cloud water content, precipitation rate and cloud optical properties affect the gas and aerosol concentration predicted by CHIMERE through aqueous chemistry, wet scavenging and photolysis.

Model Geometry
As shown in Figure 1, the model was configured with three one-way nested domains. The outer domain (D1) has 36 km resolution, domain 2 (D2) 12 km, and inner domain (D3) is cloud resolving at 4 km resolution. D1 extends approximately from latitudes 15° N to 55° N and from longitudes 15° W to 32° E with 150 × 150 grid points. D2 and D3 are centered on the Benelux Union and include 121 × 121 and 166 × 166 nodes, respectively. The reason to have a very wide-ranging domain as D1 is motivated by the need to have an accurate transport of mineral dust because, as will see in Section 4.1, our simulations on D3 are affected by long-range transport of Saharan aerosol dust for a relatively long period. In this study, we used 33 vertical levels extending from surface up to 50 hPa for WRF, and 20 levels up to 200 hPa for CHIMERE. Above CHIMERE top, aerosol concentration exchange between the models are assumed to decrease with an exponential profile. The coupling time (exchange frequency) between models is resolution dependent, it is fixed at 10 min for D1, 6 min for D2, and 5 min for the high resolution domain D3. The coupling frequency increases with the horizontal resolution in order to capture most of the meteorological variability [51].

Simulation Period
The period simulated in this study extends from 1 to 30 May 2008. This period is splitted into two different periods: 1.
From 1 to 13 May 2008: this is the spin-up period and the simulation results will not be analyzed.

2.
From 14 to 30 May 2008: this is the period where the simulation results will be compared to the available measurement. During this analyzed period, different meteorological periods were observed. This leads to different behaviour of the transport, chemistry and, thus, the potential indirect effects of aerosol on clouds. To better understand the following results, this period is divided in four sub-periods as: For these two weeks, a series of 24 h simulations was performed each day starting at 00 UTC. Each simulation is restarted from a 12 h of meteorological pre-forecast starting at 12 UTC the previous day in order to spin-up the meteorology. The choice to re-initialize the WRF model every 24 h instead of performing a continuous run is motivated by the results obtained from preliminary sensitivity tests conducted with WRF offline (no feedback with CHIMERE). These tests have shown that the re-initialization improves the simulation of cloud fields in the cloud resolving domain with respect to a continuous run.

Parametrizations Used
Parameterizations for meteorological and chemical processes used in this work are listed in Table 1. Main parameterizations used in WRF include the Thompson and Eidhammer [40] aerosol aware cloud microphysics scheme modified as described in Section 2, RRTMG for visible and planetary radiation [39], YSU for planetary boundary layer [52], Noah land surface model [53], and Grell-Freitas [54] for cumulus convection. No cumulus parameterization is used in the cloud resolving domain.
Gas-phase chemical evolution in CHIMERE is simulated using the MELCHIOR2 mechanism [55], a reduced version of MELCHIOR1 [56], including 40 gaseous species and 120 photochemical reactions. Photolytic rates are calculated with the Fast-JX scheme [57,58] as described by Mailler et al. [59]. Aerosol particles in CHIMERE are parameterized with the sectional aerosol model described by Bessagnet et al. [60]. In this study, 10 dimensional bins are used. The aerosol species treated by CHIMERE are sulfate (SO 4 ), nitrate (NO 3 ), ammonium (NH 4 ), black carbon (BC), primary unspeciated aerosol, organic matter (OM) constituted by primary and secondary organic aerosols (POA and SOA), aerosol water, mineral dust and sea salt. SOAs are formed by anthropogenic and biogenic sources as described in Bessagnet et al. [61]. Aerosols are assumed to be in thermodynamic equilibrium calculated with a lookup table as in Menut et al. [62]. Dry deposition velocities of trace gases and aerosols are determined according to Wesely [63] and Zhang et al. [64], respectively. The treatment of in and below clouds wet deposition of gases and aerosols follows Menut et al. [62].
Anthropogenic emissions are from the Hemispheric Transport of Air Pollution (HTAPv2) [65] global inventory as described in Mailler et al. [36]. Biogenic, dust and sea salt emissions are calculated online within CHIMERE. Biogenic emissions are estimated with the model of Emissions of Gases and Aerosols from Nature (MEGAN) [66]. Mineral dust and sea salt emissions are parameterized following Menut et al. [67] and Monahan [68], respectively. Initial conditions for WRF were provided by 6-hourly operational NCEP analyses at a resolution of 1 • × 1 • . Meteorological boundary conditions for Domain 1 were taken from NCEP data. On the boundaries of D1, CHIMERE is forced with the climatology of gases and aerosols taken from output of global model LMDZ-INCA [69], while the GOCART model [70] is used for dust boundary conditions. For the inner domains D2 and D3, meteorological and chemical boundary conditions were taken from D1 and D2 simulations, respectively. For all three domains, the chemical initial conditions are restarted from previous run. To improve the meteorological prediction, the spectral nudging toward NCEP analysis of wind, temperature and geopotential has been applied in D1 above 850 hPa with a relaxation time of about 1 h.

Model Evaluation
WRF-CHIMERE model including the aerosol indirect feedback is evaluated in cloud resolving domain by using aerosol measurements issued in May 2008 during the Intensive Cloud Aerosol Measurement Campaign (IMPACT) campaign in the frame of European Integrated project on Aerosol Cloud Climate and Air quality interaction (EUCAARI) project [35]. Model results are also evaluated through cloud products of MODIS satellite.
EUCAARI dataset is a state-of-the-art collection of aerosol measurements available in Europe, it is a unique and complete dataset available for model evaluation of ACI and in the past was already used in many modelling studies (e.g., [71][72][73][74][75][76]).

Meterological Conditions
From a meteorological point of view, May 2008 is divided in two periods. These two periods led to different conditions of aerosol load and composition. During the first two weeks, Europe was dominated by an anticyclonic blocking, whereas the second part of the month was characterized by the passage of frontal systems over Central Europe. An overview of the synoptic situation is provided by Hamburger et al. [77].
As shown in Figure S1, on 14 May a high pressure system was centered over Iceland, maintaining stable weather over Benelux where the 850 hPa wind circulation was from South-East because of a cyclone located over the Atlantic. Later the same day, the anticyclone started to decay and the day after the Atlantic cyclone reached continental Europe through the Bay of Biscay ( Figure S2).
On 16 May a depression moved from Scandinavia toward Central Europe and it merged with the Atlantic depression in the first hours of 17 May ( Figure S3). Several frontal systems crossed the Benelux Union, resulting in an advection of cold and clean Arctic air through the North Sea towards Central Europe. Starting from 18 May the onset of a North-Esterly anticyclonic field and the entry of Atlantic unstable air on the Mediterranean basin favoured the formation of large depressions over north Italy ( Figure S4). The cyclone located over north Italy advected warm and moist air towards Central Europe. At the same time, the Scandinavian depression advected cold air over Central Europe through the North Sea. The convergence line of these different air masses caused unstable weather on the Central Europe until 22 May. The days between 16 and 22 May are referred to as a "scavenging background period" [78].
On 21 May, a cut-off broked-off from Iceland cyclone and reached the Biscay Bay on 24 May forming a deep depression until 30 May ( Figure S5). During this period, the resulting mean circulation on the Benelux Union was characterized by a warm and moist Southerly flow that allowed a long-range transport of dust from the Sahara desert [74].

Ground Based and Aircraft Observations
Aerosol mass measurements were performed at Cabauw Experimental Site for Atmospheric Research (CESAR) observatory. It is a tower located 50 km south of Amsterdam ( Figure 2). Aerosol mass concentration of SO 4 , NO 3 , NH 4 , and OM is measured with an aerosol mass spectrometry (AMS) located at 60 m on CESAR tower [78]. In addition, surface measurements of daily inorganic mass concentration, PM 2.5 , and PM 10 are also available with the EMEP and AIRBASE European networks.
through the North Sea. The convergence line of these different air masses caused unstable weather on the Central Europe until 22 May. The days between 16 and 22 May are referred to as a "scavenging background period" [78].
On 21 May, a cut-off broked-off from Iceland cyclone and reached the Biscay Bay on 24 May forming a deep depression until 30 May ( Figure S5). During this period, the resulting mean circulation on the Benelux Union was characterized by a warm and moist Southerly flow that allowed a long-range transport of dust from the Sahara desert [74].

Ground Based and Aircraft Observations
Aerosol mass measurements were performed at Cabauw Experimental Site for Atmospheric Research (CESAR) observatory. It is a tower located 50 km south of Amsterdam ( Figure 2). Aerosol mass concentration of SO4, NO3, NH4, and OM is measured with an aerosol mass spectrometry (AMS) located at 60 m on CESAR tower [78]. In addition, surface measurements of daily inorganic mass concentration, PM2.5, and PM10 are also available with the EMEP and AIRBASE European networks.
Vertical profiles of aerosol mass concentration, aerosol size distribution, and CCN were collected aboard the French ATR-42 aircraft. In this study, we used the flights performed from 15 to 30 May 2008 within D3. Figure 2 shows the flight tracks; date, number, takeoff and landing time are reported in Table 2. A detailed description of the aircraft campaign and the ATR-42 instrumentation is provided by Crumeyrolle et al. [79]. Aerosol mass concentration of several aerosol species was measured with an AMS, a cloud condensation nuclei counter (CCNC) was used to measure the CCN concentration. Aerosol size distribution of particle in the range of 0.02-0.5 μm was measured with a scanning mobility particle sizer (SMPS), while number size distribution of aerosol larger than 0.1 μm was measured with a passive cavity aerosol spectrometer probe (PCASP).  Vertical profiles of aerosol mass concentration, aerosol size distribution, and CCN were collected aboard the French ATR-42 aircraft. In this study, we used the flights performed from 15 to 30 May 2008 within D3. Figure 2 shows the flight tracks; date, number, takeoff and landing time are reported in Table 2. A detailed description of the aircraft campaign and the ATR-42 instrumentation is provided by Crumeyrolle et al. [79]. Aerosol mass concentration of several aerosol species was measured with an AMS, a cloud condensation nuclei counter (CCNC) was used to measure the CCN concentration. Aerosol size distribution of particle in the range of 0.02-0.5 µm was measured with a scanning mobility

Cloud Optical Properties
Model skill in reproducing cloud optical properties was assessed by comparing WRF-CHIMERE output with MODIS L2 cloud products. Model output was sampled at the time of MODIS overpass on the high resolution domain, D3. MODIS and modelled data were regridded and averaged on a same regular grid with horizontal resolution of 5 km. The "Cloud_Phase_Optical_Properties" MODIS product has been used to distinguish liquid and ice cloud [80].

Metric for Model Evaluation
All ground-based and aircraft measurements described in the previous section have been used to evaluate WRF-CHIMERE performances within the high resolution domain D3. The model skill in reproducing the observations was quantified through the use of some statistical indices as the trimean (TM), the interquantile range (IQR), correlation coefficient (r), mean absolute bias (MB), normalized mean bias (NMB), and normalized mean gross error (NMGE). The indices are defined in Appendix A. The use of TM and IQR instead of the classic mean and standard deviation prevents the model evaluation from being affected by outliers [81]. IQR quantifies the spread of the data around the central tendency given by TM.

Aerosol Mass Concentration
In this subsection we assess the WRF-CHIMERE skill in reproducing the aerosol mass concentration observed at surface and aboard the ATR42. Figure 3 shows the comparison between the observed and predicted time series of daily TM of SO 4 , NO 3 , NH 4 and OM at Cabauw Tower, with shaded area and red bar denoting the 25th and 75th percentiles. Daily TMs and interquartile ranges are calculated starting from hourly modelled and observed and modelled values. Statistical indices obtained from the comparison are listed in Table 3.
Atmosphere 2019, 10, x FOR PEER REVIEW 9 of 31 are calculated starting from hourly modelled and observed and modelled values. Statistical indices obtained from the comparison are listed in Table 3.  In general, the model reproduces the day-to-day variations of observations, especially for inorganic aerosols. WRF-CHIMERE captures the marked decrease of the concentrations during the scavenging period, the following recovery starting from 20 May, and the sharp increase of aerosol pollution at the end of month. Model performance is poor on 15-16 May, all three modelled aerosol inorganic species exhibit a peaks 4-5 times larger than observations. This unrealistic modelled peak is not representative of the chemistry and transport ability of the model. For all other days, the model is able to reproduce the observations. For this peak in particular, the fact to have an over-estimation for the three modelled inorganic species and only for one day, shows that the problem is not coming from emissions or chemistry, but more probably to a badly reproduced local wind speed or direction, advecting a polluted air mass just over the tower.  Table 3. Trimean values (TM), interquantile range (IQR), correlation coefficient (r), mean bias (MB), normalized mean bias (NMB), and normalized mean gross error (NMGE) from the comparison between observed and modelled aerosol sulfate (SO 4 ), aerosol nitrate (NO 3 ), aerosol ammonium (NH 4 ), and organic matter (OM). In general, the model reproduces the day-to-day variations of observations, especially for inorganic aerosols. WRF-CHIMERE captures the marked decrease of the concentrations during the scavenging period, the following recovery starting from 20 May, and the sharp increase of aerosol pollution at the end of month. Model performance is poor on 15-16 May, all three modelled aerosol inorganic species exhibit a peaks 4-5 times larger than observations. This unrealistic modelled peak is not representative of the chemistry and transport ability of the model. For all other days, the model is able to reproduce the observations. For this peak in particular, the fact to have an over-estimation for the three modelled inorganic species and only for one day, shows that the problem is not coming from emissions or chemistry, but more probably to a badly reproduced local wind speed or direction, advecting a polluted air mass just over the tower.

Species
However, correlation between modelled and measured TM is larger than 0.5 for all aerosol species analyzed, it is 0.55, 0.54, 0.52 and 0.64 for SO 4 , NO 3 , NH 4 and OM, respectively. Despite the overestimation of 15-16 May, NMB of WRF-CHIMERE TM is relatively small for aerosol sulfate (+14%) and ammonium (+11%). Instead, particulate nitrate is overestimated by 48%. TM of OM is underestimated over all simulated periods, MB is −1.7 µg/m 3 (or −37%). The analysis of IQR reveals that model spread around central tendency of the distribution is about 1.6 times larger than observed for inorganic aerosols. By contrast, it is 1.7 times lower for OM.
The model tendency to overestimate surface aerosol nitrate is also confirmed by comparison of WRF-CHIMERE simulations with observations of daily inorganic aerosol mass concentration observed at EMEP sites. As shown in Table 4, daily NO 3 is highly biased by 0.6 µg/m 3 (+247%). As at Cesar tower, modelled SO 4 shows a relative small bias of +14%. Daily aerosol ammonium is overpredicted by 0.2 µg/m 3 (+44%). SO 4 , NO 3 and NH 4 are reproduced with a correlation of 0.48, 0.70, and 0.56, respectively. These results are statistically comparable with those obtained by large-scale simulation in other modelling studies over Europe. For example, in multi-year simulations of Lecoeur and Seigneur [82] correlation coefficient is 0.57, 0.42, and 0.58 for SO 4 , NO 3 and NH 4 , respectively. Balzarini et al. [83] reported correlations of 0.48, 0.60, and 0.56 for the same variables. For the same studied case, Tuccella et al. [76] found that daily SO 4 , NO 3 and NH 4 are simulated with a correlation of 0.66, 0.74 and 0.82, respectively. All these results are using modeling systems very different compared to the one used in this study, and it is difficult to associate the model set-up, the meteorology or the emissions to the statistical results found. The main message is that the model presented in this study is able to have realistic results and with the same order of accuracy than other existing online modeling systems.  Table 4 also reports the statistical indices related to daily PM 2.5 and PM 10 at EMEP stations within D3. WRF-CHIMERE reproduces the day-to-day variability of PM 2.5 and PM 10 with a correlation of 0.54 and 0.51, and mean positive bias of 2.4 and 5.6 µg/m 3 , respectively. Time series (not shown) show that PM 2.5 is systematically overestimated during the whole period, while for PM 10 the bias is essentially due to the long-range transport of dust. It should be noted that, as shown by Mailler et al. [36] in a continental scale evaluation of CHIMERE, the current version of the model tends to overestimate the PM 2.5 in both winter and summer. Mailler et al. [36] associated the bias to ammonium nitrate overprediction, and this is coherent with the results discussed above.
The analysis of vertical profiles has been diversified for the planetary boundary layer (PBL) and free troposphere (FT). According to Crumeyrolle et al. [80], PBL height was lower 1600 during the EUCAARI campaign. Following Tuccella et al. [76], we have considered as PBL concentrations the modelled data below 1600 m and as FT values those above this altitude. Figure 4 shows the comparison between the predicted and observed vertical profiles of SO 4 , NO 3 , NH 4  During P1, WRF-CHIMERE captures the profile of observed mass concentration in the marine PBL, except for SO4 that is underestimated by more than a factor of 2 close to the surface. Moreover, during P1, the model exhibits a considerable underestimation above the Planetary Boundary Layer (PBL). Measurements above the PBL during P1 reported in this study are not representative of clean marine air, because they were mainly carried out close to the Dutch coast. The bias above the PBL is most probably due to errors in the wind simulation since the simulated plume is localized further to the South compared to observations (not shown).
During P2a and P2b, the model tends to underestimate the mass concentration of observed SO4, NO3, NH4 and OM especially close to the surface. Moreover, WRF-CHIMERE does not reproduce the During P1, WRF-CHIMERE captures the profile of observed mass concentration in the marine PBL, except for SO 4 that is underestimated by more than a factor of 2 close to the surface. Moreover, during P1, the model exhibits a considerable underestimation above the Planetary Boundary Layer (PBL). Measurements above the PBL during P1 reported in this study are not representative of clean marine air, because they were mainly carried out close to the Dutch coast. The bias above the PBL is most probably due to errors in the wind simulation since the simulated plume is localized further to the South compared to observations (not shown).
During P2a and P2b, the model tends to underestimate the mass concentration of observed SO 4 , NO 3 , NH 4 and OM especially close to the surface. Moreover, WRF-CHIMERE does not reproduce the observed variability especially during the period 2a. The best agreement between observations and model results was found during P3. During this period, modelled values are within the observed range except for the layers close to the surface. WRF-CHIMERE captures several features of observed vertical shape. Both observed and modelled aerosols increase in the PBL, reach the maximum at an altitude of about 500 m, and show a homogeneous profile above the PBL. Figure 5 shows the comparison between the observed and predicted vertical profile of PM 2.5 . During P1, PM 2.5 is within the observed range in the PBL but is underestimated by a factor 2-3 in the Free Troposphere (FT). The comparison with the flights carried out during the scavenging period (P2a and P2b) shows that the modelled PM 2.5 in the PBL is at the lower end of the observed distribution, but the agreement is satisfying in the FT. Finally, during P3, WRF-CHIMERE underestimates the PM 2.5 mass concentration less than a factor 2 close to the surface, while PM 2.5 is overestimated on the rest of vertical profile. The positive bias in the FT is most likely attributable to the fine fraction of advected desert dust during P3. The overestimation of the dust fine fraction is consistent with the results reported by Menut et al. [84].   The black and red bars represent the observed and predicted 25th and 75th percentiles, respectively. Table 5 summarizes the statistical indices related to the comparison between predicted and observed aerosol mass concentration over the whole campaign in the PBL and FT. WRF-CHIMERE underestimates the TM of SO 4 on average by about 1.1 µg/m 3 (or −50%) in both PBL and FT. By contrast, NO 3 TM is overestimated by 0.5 (+26%) and 0.9 µg/m 3 (+450%) in the PBL and FT, respectively. However, looking at Figure 4, we underline that this bias is essentially due to P2a period and therefore does not reflect a systematic error of the model in aerosol nitrate simulation. Particulate ammonium is simulated within a factor 2 with respect to the observations, the bias is −39% and −30% in the PBL and FT, respectively. OM is systematically underestimated by a factor 2-3 in both PBL and FT. Average PM 2.5 is underpredicted by the model by −6 µg/m 3 (−32%) in the PBL, but shows a negligible error (+7%) in the FT. The analysis of IQR reveals the general tendency of the model to underestimate the variability of aerosol mass concentration found in the observations, except for NO 3 and NH 4 in the PBL and FT respectively. The underestimation of the aerosol mass concentration variability is a common feature found in other modeling systems (e.g., Tuccella et al. [76]). Part of the model underestimation during P2a and P2b could be related to an excess of wet scavenging in The underestimation of the aerosol mass concentration variability is a common feature found in other modeling systems (e.g., Tuccella et al. [76]). Part of the model underestimation during P2a and P2b could be related to an excess of wet scavenging in the model. Another critical factor for a correct simulation of bulk aerosol mass is a proper simulation of the aerosol size distribution. This will be discussed in the Section 4.5.2.

Aerosol Size Distribution and CCN
The first step for reliable simulations of ACI is a satisfactory prediction of the aerosol size number distribution together with a proper representation of aerosol mass and composition. These three elements control the number of aerosol particles activated as cloud droplet. Figure 6 shows the TMs of aerosol size distribution simulated and observed aboard the ATR-42 in all three periods within the PBL and FT. WRF-CHIMERE exhibits a larger number of aerosol particles in the Aitken mode with respect to the observations in both PBL and FT during P1, P2, and P3. CHIMERE particles relative to the bins of Aitken mode (aerosol particle with diameter between 0.01 and 0.1 µm) are overestimated from a few hundred up to 2-3 orders of magnitude.
There is also a tendency to underestimate the larger particles of the accumulation mode. These results suggest that the overestimation of ultrafine particles is likely related to an excess of particle formation from nucleation. Moreover, the lack of sufficient particles in the accumulation mode reduces the coagulation of ultrafine aerosols on larger existing particles, so that the loss of smaller particles by coagulation is not efficient enough. Another probable source of error is the size distribution at emission: a profile is applied to the total mass of primary aerosol emissions into the model bins (fixed for each type of source and the emitted species). CHIMERE does not include a prognostic treatment of particle number, which is diagnosed from the mass, density and particle diameter. It should be noted that ACI is unlikely affected by ultrafine particles, the most favoured particles to act as CCN are those larger than 100 nm. Therefore, it is interesting to explore how the model reproduces these particles referred to as condensation nuclei (CN) in this study. There is also a tendency to underestimate the larger particles of the accumulation mode. These results suggest that the overestimation of ultrafine particles is likely related to an excess of particle formation from nucleation. Moreover, the lack of sufficient particles in the accumulation mode reduces the coagulation of ultrafine aerosols on larger existing particles, so that the loss of smaller particles by coagulation is not efficient enough. Another probable source of error is the size distribution at emission: a profile is applied to the total mass of primary aerosol emissions into the model bins (fixed for each type of source and the emitted species). CHIMERE does not include a prognostic treatment of particle number, which is diagnosed from the mass, density and particle diameter. It should be noted that ACI is unlikely affected by ultrafine particles, the most favoured particles to act as CCN are those larger than 100 nm. Therefore, it is interesting to explore how the model reproduces these particles referred to as condensation nuclei (CN) in this study. Table 6 reports the observed and predicted TMs of CN within the PBL and FT. CN are underpredicted in the PBL from a factor 2.8 during P1 to negligible bias during P2a and P2b. By contrast, WRF-CHIMERE predicts a larger number of CN (about a factor 2) in the PBL during the large transport days when a significant contribution to CN is given by the fine fraction of the dust advected from the desert. In the FT, the model exhibits a systematic positive bias ranging from a factor 1.3 up to 3.3. Looking at Figure 6, it should be observed that the bias in simulating CN is not uniformly distributed with the aerosol diameter but in both case of low and high bias, the model in general overestimates the CN with diameter between 0.1 and 0.2 μm and underestimates the particles larger than 0.2 μm. This bias could have an important impact on the CCN calculation leading to larger (lower) concentrations at high (low) levels of supersaturation.
Vertical profiles of CCN at 0.2% of supersaturation are shown in Figure 7 and their TMs for each period are reported in Table 6. In the PBL, CCN concentration is underestimated during wet scavenging days and overestimated during P1 and P3, but the bias is always lower than a factor 1.5. Concerning the FT, although CN are high biased in each period, modelled CCN are lower relative to the observations by a factor 3.4 during P2b, 1.6 during P2a and P3, respectively. CCN are  Table 6 reports the observed and predicted TMs of CN within the PBL and FT. CN are underpredicted in the PBL from a factor 2.8 during P1 to negligible bias during P2a and P2b. By contrast, WRF-CHIMERE predicts a larger number of CN (about a factor 2) in the PBL during the large transport days when a significant contribution to CN is given by the fine fraction of the dust advected from the desert. In the FT, the model exhibits a systematic positive bias ranging from a factor 1.3 up to 3.3. Looking at Figure 6, it should be observed that the bias in simulating CN is not uniformly distributed with the aerosol diameter but in both case of low and high bias, the model in general overestimates the CN with diameter between 0.1 and 0.2 µm and underestimates the particles larger than 0.2 µm. This bias could have an important impact on the CCN calculation leading to larger (lower) concentrations at high (low) levels of supersaturation.
Vertical profiles of CCN at 0.2% of supersaturation are shown in Figure 7 and their TMs for each period are reported in Table 6. In the PBL, CCN concentration is underestimated during wet scavenging days and overestimated during P1 and P3, but the bias is always lower than a factor 1.5. Concerning the FT, although CN are high biased in each period, modelled CCN are lower relative to the observations by a factor 3.4 during P2b, 1.6 during P2a and P3, respectively. CCN are overestimated in the FT only during P1. The results found here are comparable with the performances of other state-of-art models. For example, Tuccella et al. [76] comparing WRF/Chem to our same data set, reported an overestimation of CCN by a factor 1.5 in the mixed layer and a small bias in the FT.  . Satellite CDNC has been calculated starting from liquid cloud droplet effective radius (Re) and liquid cloud optical depth (LCOD) following the method used by Georg and Wood [85]. WRF-CHIMERE captures qualitatively several features of the spatial distribution of observed CDNC. The gradient between North Sea and land is reproduced by model, as well as the distribution of the maxima values of CDNC, with the exception of some regions where the model shows an important bias, for example South-East of England coast, in the middle of the border between the Netherlands and Germany, on the border between France and Belgium. As shown in Table 7, domain average of calculated CDNC is underestimated by 40%.
As it is evident from Figure 8a,b, both model and observations show a fine scale variability that makes a point-by-point comparison between WRF-CHIMERE and MODIS difficult. Modelled and observed clouds exhibit small differences in location and timing of formation. As a consequence, the distribution function (DF) is a useful tool that provides an overall view of WRF-CHIMERE performance in reproducing CDNC and at the same time is independent by spatial and temporal In addition to the size distribution, modelled CCN bias is associated with the critical diameter (CD) and aerosol composition. CD is defined as the cut-off diameter above which all particles act as CCN at a given level of supersaturation. Observed CD has been calculated from aerosol size number distribution and CCN measurements. It is given by the diameter at which the integrated size distribution equals the CCN concentration. TMs of observed and predicted CD are reported in Table 6 and are overlaid to the size distributions shown in Figure 6. To explain the link between aerosol size distribution, CD and composition in reproducing CCN, we focus as an example on the results obtained in the FT during P3. We choose this case because, as shown in Figure 6, the predicted size distribution is very similar to the observations for diameters larger than about 0.15 µm, therefore the CCN bias has a negligible dependence on the size distribution. In the case examined, CD is larger than the observed one by a factor 1.5, explaining the negative bias of CCN despite the reasonable simulation of the size distribution (for diameters larger than 0.15 µm). According to Kohler theory, CD depends on the aerosol hygroscopicity that in turn is a function of composition. Given a supersaturation level, the larger the hygroscopicity, the lower the CD. As a consequence, underestimation of hygroscopicity in the FT during P3 is most likely due to an excess of dust in the fine fraction of modelled PM 2.5 (see Section 4.5.1) combined with the missing of sufficient inorganic (sulphate, ammonium, nitrate) mass concentration (about 50% less compared to AMS data) that lowers the aerosol solubility. A measure of hygroscopicity is given by CCN efficiency (CCNE) defined as the ratio between CCN and CN. The observed and modelled values are shown in Table 6. As expected, CCNE simulated by WRF-CHIMERE for the case examined is about 5 times smaller than the observations. In the PBL, CCNE is simulated within a factor 2 but the error is larger the FT.   In general, WRF-CHIMERE overestimates the occurrence frequency (OF) of CDNC less than 100 #/cm 3 and underpredicts the largest values. Specifically, the cumulated OF of CDNC representative of "clean atmospheres" (less than 100 #/cm 3 ) are overestimated by a factor of 2. By contrast, modelled CDNC typical of "polluted conditions" (100-600 #/cm 3 ) are about 1.5 times smaller relative to the MODIS data, while values representative of "very polluted conditions" (larger than 600 #/cm 3 ) constitute the 0.5% and 1.5% of observed and modelled sample, respectively.  As it is evident from Figure 8a,b, both model and observations show a fine scale variability that makes a point-by-point comparison between WRF-CHIMERE and MODIS difficult. Modelled and observed clouds exhibit small differences in location and timing of formation. As a consequence, the distribution function (DF) is a useful tool that provides an overall view of WRF-CHIMERE performance in reproducing CDNC and at the same time is independent by spatial and temporal phase errors [86]. Observed and modelled DFs of CDNC are shown in Figure 8c.

Cloud Optical Properties
In general, WRF-CHIMERE overestimates the occurrence frequency (OF) of CDNC less than 100 #/cm 3 and underpredicts the largest values. Specifically, the cumulated OF of CDNC representative of "clean atmospheres" (less than 100 #/cm 3 ) are overestimated by a factor of 2. By contrast, modelled CDNC typical of "polluted conditions" (100-600 #/cm 3 ) are about 1.5 times smaller relative to the MODIS data, while values representative of "very polluted conditions" (larger than 600 #/cm 3 ) constitute the 0.5% and 1.5% of observed and modelled sample, respectively.
The origin of the bias of modelled CDNC is related to the uncertainties in aerosol size distribution and hygroscopicity prediction. As already mentioned in the Section 4.5, CCN are underestimated at low supersaturations and high biased for large supersaturations, because modelled CN are larger than observations at the smallest diameters of accumulation mode (0.1-0.2 µm) but are low biased for larger sizes. Therefore, given that the range of supersaturation of our simulations is relatively low, the main uncertainties in prediction of CDNC is due to the general underestimation of CN at larger diameter of aerosol size distribution. Another factor that explains CDNC bias is the uncertainty related to hygroscopitcity that in turn affects the CD. As discussed in the previous Section, modelled CD is larger with respect to the observed one because CCNE (which is a measure of aerosol hygroscopitcity) predicted by WRF-CHIMERE is low relative to ATR-42 data. Larger CD implies a lower amount of CCN that may be activated as cloud droplet. Figure 9a,b shows the TMs calculated over the whole period of liquid cloud droplet effective radius (R e ). As for CDNC, WRF-CHIMERE is able to reproduce qualitatively the features of observed Re pattern. The gradient between North Sea and land and the distribution of the maxima values of Re are reproduced by the model, with the exception of the regions where a bias in CDNC. In these regions, CDNC was noticed and Re are anticorrelated. As shown in Table 7, modelled and measured TM and IQR of R e show very similar values, meaning that WRF-CHIMERE is able to capture the central tendency of measured R e .
As displayed in Figure 9c, the maximum of DF of modelled R e has the peak shifted toward larger values with respect to satellite data. DFs of MODIS and WRF-CHIMERE exhibit the maximum around 9-10 and 7-8 µm, respectively. Cumulated OF around the observed maximum (6-9 µm) represents about the 50% of MODIS sample, while in the same range falls the 34% of WRF-CHIMERE data. The shift of the peak toward larger droplet found in WRF-CHIMERE is a consequence of low CDNC in the model. Figure 10a,b report the TMs of retrieved and predicted LCOD, respectively. As for previous variables, WRF-CHIMERE explains the features of LCOD observed by MODIS. It reproduces the lowest values of the pattern observed in the North of the domain, and the band of maximum values extending from northern France up to northwest Germany crossing Belgium and The Netherland.
WRF-CHIMERE overestimates observed LCOD along this band, but the model tends to slightly underpredict LCOD elsewhere. LCOD averaged over the domain is reported in Table 7. It is underestimated by a factor 1.3 (−19%) but exhibits larger variability (IQR) than the observations. As displayed in Figure 9c, the maximum of DF of modelled Re has the peak shifted toward larger values with respect to satellite data. DFs of MODIS and WRF-CHIMERE exhibit the maximum around 9-10 and 7-8 μm, respectively. Cumulated OF around the observed maximum (6-9 μm) represents about the 50% of MODIS sample, while in the same range falls the 34% of WRF-CHIMERE data. The shift of the peak toward larger droplet found in WRF-CHIMERE is a consequence of low CDNC in the model. Figure 10a,b report the TMs of retrieved and predicted LCOD, respectively. As for previous variables, WRF-CHIMERE explains the features of LCOD observed by MODIS. It reproduces the lowest values of the pattern observed in the North of the domain, and the band of maximum values extending from northern France up to northwest Germany crossing Belgium and The Netherland. WRF-CHIMERE overestimates observed LCOD along this band, but the model tends to slightly underpredict LCOD elsewhere. LCOD averaged over the domain is reported in Table 7. It is underestimated by a factor 1.3 (−19%) but exhibits larger variability (IQR) than the observations. As we would expect, the bias in simulation of cloud optical properties depend on the particular aerosol mechanism adopted in the models. For example, Tuccella et al. [76], using WRF/Chem with a modal aerosol model on the same region and period of our study, found an opposite behavior with respect to our results. Indeed, Tuccella et al. [76] reported an overestimation of CCN that produces smaller Re and larger LCOD relative to the observations. However, it should be noted that the magnitude of bias found here is comparable to WRF/Chem performances shown by Tuccella et al. [76], for example LCOD is overestimated by a few percent up to 48%. As described in Section 2, WRF-CHIMERE includes also the glaciation indirect effect through the processes of heterogeneous and homogeneous ice nucleation starting from aerosol predicted by chemical modules. Unfortunately, no available in situ measurements allow the evaluation of the ice particle concentration. The simulated ice cloud optical depth (ICOD) was compared to that observed by MODIS, including the precipitating hydrometeors and suspended ice cloud. As we would expect, the bias in simulation of cloud optical properties depend on the particular aerosol mechanism adopted in the models. For example, Tuccella et al. [76], using WRF/Chem with a modal aerosol model on the same region and period of our study, found an opposite behavior with respect to our results. Indeed, Tuccella et al. [76] reported an overestimation of CCN that produces smaller Re and larger LCOD relative to the observations. However, it should be noted that the magnitude of bias found here is comparable to WRF/Chem performances shown by Tuccella et al. [76], for example LCOD is overestimated by a few percent up to 48%.
As described in Section 2, WRF-CHIMERE includes also the glaciation indirect effect through the processes of heterogeneous and homogeneous ice nucleation starting from aerosol predicted by chemical modules. Unfortunately, no available in situ measurements allow the evaluation of the ice particle concentration. The simulated ice cloud optical depth (ICOD) was compared to that observed by MODIS, including the precipitating hydrometeors and suspended ice cloud. Figure 10c,d shows the TMs of observed and simulated ICOD. WRF-CHIMERE captures qualitatively the maximum values displayed in the observed pattern in a band between 52-53 • N and on the border between Belgium and Germany. WRF-CHIMERE tends to underpredict the minimum values of ICOD observed over the sea and Northwestern Germany. Simulated ICOD TM average is about half that of observations, but IQRs are very similar meaning larger dispersion of modelled values with respect to MODIS. Figure 11b reports the comparison between observed and modelled ICOD DFs. The comparison suggests that WRF-CHIMERE has about 15% more of ICOD less than 2. DFs reach maximum for values between 3-20 and 20-40 for observed and predicted sample, respectively. WRF-CHIMERE underestimate the maximum of observed OF by about 15-20%. As described in Section 2, WRF-CHIMERE includes also the glaciation indirect effect through the processes of heterogeneous and homogeneous ice nucleation starting from aerosol predicted by chemical modules. Unfortunately, no available in situ measurements allow the evaluation of the ice particle concentration. The simulated ice cloud optical depth (ICOD) was compared to that observed by MODIS, including the precipitating hydrometeors and suspended ice cloud. Figure 10c,d shows the TMs of observed and simulated ICOD. WRF-CHIMERE captures qualitatively the maximum values displayed in the observed pattern in a band between 52-53° N and on the border between Belgium and Germany. WRF-CHIMERE tends to underpredict the minimum values of ICOD observed over the sea and Northwestern Germany. Simulated ICOD TM average is about half that of observations, but IQRs are very similar meaning larger dispersion of modelled values with respect to MODIS. Figure 11b reports the comparison between observed and modelled ICOD DFs. The comparison suggests that WRF-CHIMERE has about 15% more of ICOD less than 2. DFs reach maximum for values between 3-20 and 20-40 for observed and predicted sample, respectively. WRF-CHIMERE underestimate the maximum of observed OF by about 15-20%.  Although in this section we have focused on the relationship among aerosol size distribution, hygroscopicity, CDNC, and COD, the uncertainties in simulating cloud optical properties may be also related to model setup used. For example, Barò et al. [87] investigated the sensitivity of WRF/Chem of ACI to two different cloud microphysics schemes, and found serious differences in CDNC and cloud water mixing ratio patterns predicted by two parameterizations that depend on geographical location, pollution level, and season. In addition to cloud microphysics representation, Otkin and Greewald [85] showed an important sensitivity of cloud optical properties to PBL schemes.

Sensitivity Test
In this section, as first example of application of WRF-CHIMERE model including ACI, we discuss the regional response to aerosol emission abatement. The studies of global warming mitigation strategies are usually conducted focusing only the reduction of emissions of aerosol species that cause the warming of the atmosphere such as the black carbon [28]. On a practical level, it is not easy to reduce the emission of single species since it is co-emitted together to other SLCPs and precursor gases that lead to secondary aerosol formation. For this reason, we did not focus on a single aerosol species, but we have reduced both primary aerosols and precursor gases emissions from all anthropogenic sources. Specifically, the method described in Section 3 was used but with anthropogenic emission fluxes halved in all three domains. Our test should be intended as a limit case of emission abatement. Aerosol direct effect is not included in the simulations. Moreover, chemical boundary conditions for D1 remain unchanged in our test. As a consequence, the effect induced by emission reduction could be underestimated. It is should be also noticed that given the short period of interest, the test is not a climatic simulation. The test is only dedicated to short-term indirect effects of aerosols on radiation, clouds, temperature and precipitations.
To investigate the change on the radiation budget due to emission reduction via indirect effects, we define the indirect radiative effect (IRE) as the average difference in the net radiation flux (down-up) at surface between sensitivity test and reference run including all current anthropogenic sources. Domain average effect of emission reduction on CCN concentration, radiation flux, and temperature is reported in Table 8. Table 8. Domain average effect of anthropogenic emission reduction. The effect on CCN and temperature is calculated as the difference between the perturbed and reference run. Indirect radiative effect (IRE) is the difference in the net radiation flux (down-up) at surface between sensitivity and perturbed simulation.

Land
Sea Land + Sea According to WRF-CHIMERE, halving of anthropogenic emissions reduces PBL average concentration of CCN at 0.2% of supersaturation by about 60%. The decrease of 60% of CCN concentration against the 50% of emission highlights the non-linear behavior of CCN formation. Figure 12 shows the average IRE calculated over the whole period (14-30 May 2008) for shortwave (SW), longwave (LW), and total radiation. Average IRE radiative effect on SW radiation is positive over whole domain and is equal to +0.89 W/m 2 . This result is consistent with our expectations. Lower CCN number concentration implies lower CDNC with larger diameter with respect to reference simulation, and thus less efficient clouds in reflecting the incoming solar radiation. Although SW IRE is positive on average, its distribution is not homogeneous on the domain but appears "noisy", locally IRE may reach values of ±12 W/m 2 . This suggests that emission abatement influences timing and location of cloud formation. Beside average IRE is different between land and sea. We notice that the amplitude of IRE is +1.56 and +0.27 W/m 2 on the land and sea respectively. On the sea, we note some regions where IRE is mainly positive.
LW IRE is also positive (+0.10 W/m 2 ) but its magnitude is less than change of SW radiation budget (about 11%). As shown in Figure 12b, LW IRE displays a marked asymmetry between land and sea. It is positive on the sea (+0.41 W/m 2 ) and negative on land (−0.24 W/m 2 ). Given the non-linearity of the aerosol-cloud-radiation system, it is not easy to explain this asymmetry. A possible explanation is related to liquid cloud top height (CTH) change that results from emission reduction. Figure 13 shows the cumulative distribution frequency of liquid CTH simulated by WRF-CHIMERE in the reference run separately for land and sea. In Figure 13 is also shown the average change in CTH for low (100 < CTH < 800 m), middle (800 < CTH < 1500 m), and high (CTH > 1500 m) liquid clouds. As we may infer from Figure 13, 80% of sea clouds is constituted by clouds with top lower than 1500, while on the land these are less than 50% of the sample. According to the model, when the anthropogenic emissions are halved, we note that top height of low and middle liquid clouds tends to increase, at the same time the top of high cloud lowers. Therefore, being the scenario on the sea dominated by low and middle cloud, the increase of their top decreases the outgoing longwave radiation. By contrast on the land where high clouds are much more frequent, the lowering of their top increase the flux of the outgoing longwave radiation. As a result, on the land SW IRE is partially compensated by LW IRE giving a net value of +1.32 W/m 2 . On the sea, SW and LW IRE have the same sign and net IRE is +0.68 W/m 2 . Domain average of IRE is +0.99 W/m 2 . LW IRE is also positive (+0.10 W/m 2 ) but its magnitude is less than change of SW radiation budget (about 11%). As shown in Figure 12b, LW IRE displays a marked asymmetry between land and sea. It is positive on the sea (+0.41 W/m 2 ) and negative on land (−0.24 W/m 2 ). Given the nonlinearity of the aerosol-cloud-radiation system, it is not easy to explain this asymmetry. A possible explanation is related to liquid cloud top height (CTH) change that results from emission reduction. Figure 13 shows the cumulative distribution frequency of liquid CTH simulated by WRF-CHIMERE in the reference run separately for land and sea. In Figure 13 is also shown the average change in CTH for low (100 < CTH < 800 m), middle (800 < CTH < 1500 m), and high (CTH > 1500 m) liquid clouds. As we may infer from Figure 13, 80% of sea clouds is constituted by clouds with top lower than 1500, while on the land these are less than 50% of the sample. According to the model, when the anthropogenic emissions are halved, we note that top height of low and middle liquid clouds tends to increase, at the same time the top of high cloud lowers. Therefore, being the scenario on the sea dominated by low and middle cloud, the increase of their top decreases the outgoing longwave radiation. By contrast on the land where high clouds are much more frequent, the lowering of their top increase the flux of the outgoing longwave radiation. As a result, on the land SW IRE is partially compensated by LW IRE giving a net value of +1.32 W/m 2 . On the sea, SW and LW IRE have the same sign and net IRE is +0.68 W/m 2 . Domain average of IRE is +0.99 W/m 2 . The effect induced by emissions reduction on temperature at 2 m (T2) and precipitation rate is mainly representative of local changes. Instantaneous changes of T2 locally may reach values of ±7 • C, while the instantaneous changes of precipitation rate are up to ±50 mm/day.
Error associated with IRE is mainly due to the uncertainties in CDNC simulation. Model evaluation has shown the model tendency to underestimate CDNC leading to larger cloud effective radius and lower cloud optical thickness with respect to the observations. As a consequence, we have to expect less incoming radiation reflected to the space by clouds (lower cloud radiative forcing) in the reference run. When CDNC is reduced in the sensitivity test, we would expect a larger warming effect compared to the one we would have in the real atmosphere, if the droplets are reduced by the same amount. Our IRE has been estimated to be overestimated by about 10%. The method used to calculate the uncertainty is reported in the Appendix B. The effect induced by emissions reduction on temperature at 2 m (T2) and precipitation rate is mainly representative of local changes. Instantaneous changes of T2 locally may reach values of ±7°C, while the instantaneous changes of precipitation rate are up to ±50 mm/day.
Error associated with IRE is mainly due to the uncertainties in CDNC simulation. Model evaluation has shown the model tendency to underestimate CDNC leading to larger cloud effective radius and lower cloud optical thickness with respect to the observations. As a consequence, we have to expect less incoming radiation reflected to the space by clouds (lower cloud radiative forcing) in the reference run. When CDNC is reduced in the sensitivity test, we would expect a larger warming effect compared to the one we would have in the real atmosphere, if the droplets are reduced by the same amount. Our IRE has been estimated to be overestimated by about 10%. The method used to calculate the uncertainty is reported in the Appendix B.

Conclusions
Aerosol indirect effects have been implemented within the online version of WRF-CHIMERE coupled model [25]. In this work, we have extended the two-way coupling to the feedback from the aerosol-cloud interaction. First, second and glaciation indirect effects have been included in the model via the cloud microphysics scheme of Thompson and Eidhammer [40]. Aerosol activation as cloud droplets, heterogeneous and homogeneous ice nucleation are calculated starting from aerosol size distribution predicted by CHIMERE model. WRF-CHIMERE with indirect effects has been evaluated in the period 14-30 May 2008 comparing cloud resolving simulations (4 km) with ground-based, aircraft measurements of aerosol mass and particle collected on the Benelux Union during the EUCAARI campaign [35] and MODIS cloud data. During the studied period, the atmospheric conditions were characterized by alternation of relative stable weather conditions, wet scavenging days, and days with long-range transport of dust from Sahara desert.
WRF-CHIMERE captures the day-to-day variations of aerosol mass concentrations at surface. Observed SO4, NO3, NH4, OM, PM2.5, and PM10 are reproduced with correlation coefficients larger than 0.50. The analysis of the bias reveals that the model tends to overestimates the observed

Conclusions
Aerosol indirect effects have been implemented within the online version of WRF-CHIMERE coupled model [25]. In this work, we have extended the two-way coupling to the feedback from the aerosol-cloud interaction. First, second and glaciation indirect effects have been included in the model via the cloud microphysics scheme of Thompson and Eidhammer [40]. Aerosol activation as cloud droplets, heterogeneous and homogeneous ice nucleation are calculated starting from aerosol size distribution predicted by CHIMERE model. WRF-CHIMERE with indirect effects has been evaluated in the period 14-30 May 2008 comparing cloud resolving simulations (4 km) with ground-based, aircraft measurements of aerosol mass and particle collected on the Benelux Union during the EUCAARI campaign [35] and MODIS cloud data. During the studied period, the atmospheric conditions were characterized by alternation of relative stable weather conditions, wet scavenging days, and days with long-range transport of dust from Sahara desert.
WRF-CHIMERE captures the day-to-day variations of aerosol mass concentrations at surface. Observed SO 4 , NO 3 , NH 4 , OM, PM 2.5 , and PM 10 are reproduced with correlation coefficients larger than 0.50. The analysis of the bias reveals that the model tends to overestimates the observed ammonium nitrate at ground stations and underestimates OM. Comparison between model and observed aircraft profiles shows that WRF-CHIMERE is able to reproduce the vertical profiles of SO 4 , NO 3 , NH 4 , OM, and PM 2.5 , but their mass concentration is underestimated in the model especially in the PBL during wet scavenging days, suggesting an excess of wet deposition. By contrast, OM is underestimated over the whole period by about a factor 2. PM 2.5 is overpredicted above the surface during the long range transport days, most likely because of an excesses of fine dust.
Model evaluation with aerosol size number distribution measured aboard the ATR-42 shows the systematic tendency of WRF-CHIMERE to overpredicted the particles of the Aitken mode in both PBL and free troposphere. Simulation of condensation nuclei (CN) (particles larger than 100 nm) does not show a systematic bias in the PBL, bias sign changes with the synoptic frame. By contrast, CN in the FT are overpredicted from a factor 1.3 up to 3.3. However, the bias in simulating CN is not uniformly distributed with the aerosol diameter but the model overestimates the CN with diameter between 0.1 and 0.2 µm and underestimates the particles larger than 0.2 µm.
Despite the uncertainties in aerosol size distribution simulation, predicted CCN vertical profiles have a good degree of similarity with the observed one. In the PBL, modelled CCN concentration is within a factor 2 with respect to the observations. Concerning the FT, although CN are high biased over the whole period, simulated CCN does not show a systematic bias, but this varies with the evolution of synoptic configuration. This behavior is explained looking at the modelled critical diameter (CD) and aerosol hygroscopicity. Taking a case where the predicted size distribution has a negligible bias with respect to the observations as for example in the FT during dust transport days, we have found that CCN efficiency (a measure of aerosol hygroscopicity defined as the ratio between CCN and CN concentrations) is 5 times smaller with respect the observations. This implies a larger CD in the model with a consequent underestimation of CCN. CCN efficiency bias is related to the uncertainties in simulating aerosol mass. In the case examined during the long range dust transport days, is most likely due to an excess of dust (hydrophobic) in the fine fraction of modelled PM 2.5 combined with the missing of sufficient inorganic salt mass (hydrophilic).
WRF-CHIMERE skill in reproducing cloud optical properties has been evaluated comparing the model results with MODIS L2 cloud products. WRF-CHIMERE reproduces several features of the average pattern of observed cloud droplet number concentration (CDNC). The cumulated observed distribution frequencies of CDNC less than 100 #/cm 3 are overestimated by a factor 2. By contrast, modelled CDNC in the range of 100-600 #/cm 3 are about 1.5 times smaller relative to the MODIS data. As for CDNC, WRF-CHIMERE is able to reproduce the average spatial distribution of observed liquid cloud effective radius (R e ) pattern. The analysis of Re distribution frequencies suggest that our model has the pick shifted toward larger values with respect to satellite data and at the same time there is a certain tendency in the model to overestimate the frequency of the smallest droplets (R e < 6 µm). Modelled average domain of liquid cloud optical depth (LCOD) is underestimated by 19% with respect to MODIS data.
The results discussed in this work indicate that the new functionality implemented in the model is able to reproduce the observations in a spatial and temporal scale dominated by high meteorological, chemical and emission variability. Compared to other similar modelling systems, results are in the same order of magnitude (although other models are never completely similar) and represent the current state of the art in regional chemistry-transport modelling. Therefore, our modeling system is suitable to be used in studies involving the feedback between aerosol and clouds. As first example of application, we have discussed the results of a sensitivity simulation where the anthropogenic emissions have been halved. The impact of emission abatement has been quantified through the indirect radiative effect (IRE) defined as the average difference in the net radiation flux (down-up) at surface between sensitivity test and reference run including all current anthropogenic sources. Average IRE on the cloud resolving domain is positive for both short and long wave radiation.
A closer inspection to the regional distribution of IRE has revealed that the change of long wave radiation budget is negative on the land and positive on the sea. This asymmetry is explainable with the change of cloud top height (CTH). The analysis of the frequency distribution of liquid CTH in the reference run has shown that cloud coverage on the sea is dominated by low (100 < CTH < 800 m) and middle (800 < CTH < 1500 m) clouds, while high clouds (CTH > 1500 m) dominate on the land. According to the model, halving of anthropogenic emissions increases the top height of low and middle liquid clouds and lowers the top of highest clouds. Therefore, on the sea the increase of low and middle CTH decreases the outgoing longwave radiation flux. By contrast on the land where high clouds are much more frequent, the lowering of their top increases the flux of the outgoing longwave radiation.
The effects induced by halving of anthropogenic emissions have been also investigated in terms of surface temperature and precipitation rate changes. Model results have shown that there are important instantaneous changes at local scale (up to ±6 • C for maximum and minimum temperatures, ±50 mm/day for precipitation). In the framework of regional air quality modelling, this could have a large impact on isolated pollutants plumes during their transport.
Although the results obtained in this work are encouraging and show that WRF-CHIMERE performances are comparable with other state-of-the-art online modeling systems, there is still much room for further evaluations and improvements. In our work WRF-CHIMERE evaluation is limited to the liquid water cloud, further verification with observations from cirrus cloud fields are recommended. Model comparison with CCN aircraft measurements and MODIS cloud data have highlighted that the uncertainties in simulating ACI in WRF-CHIMERE model manly arises from error in prediction of aerosol size number distribution and composition. For this reason, ACI simulation could improve updating aerosol model with the recently introduced scheme of Couvidat et al. [88]. Finally, ACI simulations with current version of WRF-CHIMERE are recommended with horizontal resolution less than 10 km because the indirect effect is considered only for grid-scale clouds. The extension of ACI to parameterized (convective) clouds is desirable and could be done through some aerosol aware convective parameterizations implemented in WRF like that of Grell and Freitas [54] or Berg et al. [89]. This development would allow simulating ACI at resolution coarser than cloud resolving scale.

Acknowledgments:
We are grateful to Antoinette A. Mensah, Suzanne Crumeyrolle and Alfons Schwarzenboeck for providing the EUCAARI measurements. We thank the WRF developers for making the model available to the scientific community. We thank also EMEP and AIRBASE for making freely available the ground-based measurements. where X is a generic vector, Z(X) is its standard score, and σ X is the standard deviation. • Mean Bias:

Conflicts of Interest
Mod i − Obs i • Normalized Mean Bias (NMB): Mod i − Obs i Obs i × 100 • Normalized Mean Gross Error (NMGE):

Appendix B
The estimation of indirect radiative effect (IRE) uncertainty may be quantified following Seinfeld and Pandis [90]. The SW radiative forcing at the top of the atmosphere exerted by liquid cloud associated with a CDNC variation is given by: where F 0 is the incoming solar radiation flux at the top of the atmosphere (400 W/m 2 ), T a the atmospheric transmission (0.76), A c the cloud fraction (about 0.3 in our cloud-resolving domain), and R is the cloud albedo given by: where τ is the cloud optical depth. In the sensitivity simulation, CDNC reported in Table 7 is reduced by about 20%. By assuming to have the same reduction in the real atmosphere, ∆F depends only on R. We may estimate R for both model and observations by using LCOT values reported in Table 7. The resulting average values are 0.68 and 0.61 for MODIS and WRF-CHIMERE, respectively. As a consequence, ∆F is +1.11 and +1.21 W/m 2 for observations and model, respectively. This mean that modelled ∆F is likely overestimated by +9%.