Ozone deposition impact assessments for forest canopies require accurate ozone flux partitioning on diurnal timescales

Abstract

Dry deposition is an important sink of tropospheric ozone that affects surface concentrations and impacts crop yields, the land carbon sink, and the terrestrial water cycle. Dry deposition pathways include plant uptake via stomata and non-stomatal removal by soils, leaf surfaces, and chemical reactions. Observational studies indicate that ozone deposition exhibits substantial temporal variability that is not reproduced by atmospheric chemistry models due to a simplified representation of vegetation uptake processes in these models. In this study, we explore the importance of stomatal and non-stomatal uptake processes in driving ozone dry deposition variability on diurnal to seasonal timescales. Specifically, we compare two land surface ozone uptake parameterizations – a commonly applied big leaf parameterization W89; and a multi-layer model (MLC-CHEM) constrained with observations – to multi-year ozone flux observations at two European measurement sites (Ispra, Italy, and Hyytiälä, Finland). We find that W89 cannot reproduce the diurnal cycle in ozone deposition due to a misrepresentation of stomatal and non-stomatal sinks at our two study sites, while MLC-CHEM accurately reproduces the different sink pathways. Evaluation of non-stomatal uptake further corroborates the previously found important roles of wet leaf uptake in the morning under humid conditions and soil uptake during warm conditions. The misrepresentation of stomatal versus non-stomatal uptake in W89 results in an overestimation of growing season cumulative ozone uptake (CUO), a metric for assessments of vegetation ozone damage, by 18 % (Ispra) and 28 % (Hyytiälä), while MLC-CHEM reproduces CUO within 7 % of the observation-inferred values. Our results indicate the need to accurately describe the partitioning of the ozone atmosphere–biosphere flux over the in-canopy stomatal and non-stomatal loss pathways to provide more confidence in atmospheric chemistry model simulations of surface ozone mixing ratios and deposition fluxes for large-scale vegetation ozone impact assessments.

1 Introduction

Ozone (${\mathrm{O}}_{\mathrm{3}}$) in the atmospheric surface layer is an air pollutant that is toxic to humans and plants. Ozone is removed by oceans, bare soil, and vegetated areas, which together are called dry deposition and account for $\pm \mathrm{15}$ %–20 % of the total tropospheric ozone sink . In vegetation canopies, the dominant deposition pathway is stomatal uptake, which typically accounts for 40 %–60 % of the total deposition to vegetation . Stomatal ozone uptake reduces carbon assimilation in vegetation , affects the terrestrial water cycle , and causes economic damage through reduced crop yield e.g.,. Besides stomatal uptake, ozone removal occurs via a range of non-stomatal removal mechanisms, such as uptake by the leaf exterior and soils, and in-canopy chemical removal involving nitrogen oxides (${\mathrm{NO}}_x$) or plant-emitted reactive carbon species. The contribution of these ozone removal processes to the total non-stomatal term is uncertain and displays temporal variability on diurnal to interannual timescales that is incompletely understood . Given that these non-stomatal removal processes act in parallel to the stomatal removal of ozone, the characterization and quantification of non-stomatal sinks is important for quantification of total and stomatal ozone uptake.

The contribution of different ozone uptake pathways cannot be routinely measured at the plant canopy level due to the various non-stomatal uptake pathways. Most studies infer stomatal conductance ($g_{\text{s}}$) from canopy-top micro-meteorological and eddy covariance observations using an inverted form of the Penman–Monteith equation e.g.,, although some studies apply alternative $g_{\text{s}}$ estimation methods based on gross primary production GPP; . In such observation-based studies, the non-stomatal ozone removal component ($g_{\text{ns}}$) is generally treated as the residual of the total uptake conductance ($g_{\text{c}}$; inferred based on the ozone dry deposition velocity) and $g_{\text{s}}$. However, sites with long-term ozone flux measurements are scarce , which limits the characterization of the seasonal to interannual temporal variability in the stomatal and non-stomatal components of ozone removal. Several campaign-based studies partitioned total canopy ozone fluxes by using ozone flux measurements along a vertical gradient to study the in-canopy flux divergence and to relate this to the vertical distribution of ozone sinks in the canopy , but these are limited to short timescales. Given the scarce availability of ozone deposition observations that span at least 1 year, and preferentially multiple years, quantifying temporal variability in stomatal and non-stomatal ozone deposition solely based on observations remains challenging.

Studies of ozone deposition (and its impacts) on regional to global scales rely on the application of atmospheric chemistry models and their dry deposition parameterizations. Many models treat deposition in a zero-dimensional manner and do not, or only implicitly, account for the variation in different in-canopy loss pathways as a function of environmental drivers and height within the canopy the big leaf approach,. Recent advances in the description of ozone deposition have been made by improving the simulation of stomatal conductance and the representation of various non-stomatal removal terms e.g., and in-canopy turbulence and radiation extinction . Additionally, some models account for vegetation ozone damage via effects on photosynthesis and stomatal conductance . Another class of models treats the canopy as a separate exchange regime with different biophysical and chemical conditions compared to the lowermost atmospheric layer and explicitly resolves in-canopy vertical gradients of ozone deposition and its driving variables by using multiple in-canopy layers e.g.,. Despite these advances in the representation of ozone deposition in atmospheric chemistry models, their application for ozone impact assessments remains a challenge. For example, the description of stomatal conductance is an important parameter for understanding year-to-year variability in impact metrics, such as cumulative uptake of ozone CUO; , but stomatal versus non-stomatal ozone flux partitioning in these models is uncertain. Additionally, spatiotemporal controls of ozone deposition pathways remain incompletely understood , in part owing to the scarcity of long-term ozone flux observations. Therefore, we here study the temporal controls on stomatal and non-stomatal ozone deposition pathways, and their implications for simulations of CUO, using two multi-year ozone deposition data sets as well as a big leaf and multi-layer parameterization of land surface ozone uptake.

Specifically, we investigate the added value of an explicit multi-layer canopy representation of ozone deposition MLC-CHEM – the Multi-Layer Canopy–CHemistry Exchange Model; compared to a commonly applied big leaf parameterization in terms of simulating ozone deposition pathways and ozone impact metrics. We first study long-term (seasonal to annual) and short-term (diurnal) temporal variability in ozone dry deposition to forest canopies at a pristine boreal site (Hyytiälä) and a pre-alpine site that frequently experiences high ozone concentrations (Ispra). We then evaluate the performance of a big leaf and a multi-layer representation of atmosphere–biosphere exchange in simulating ozone dry deposition pathways and their temporal variability. Subsequently, we characterize the relationship of non-stomatal conductance as a function of environmental drivers. Last, we aim to demonstrate how representations of the drivers of long- and short-term variability in ozone stomatal and non-stomatal removal in those different land surface parameterizations affect simulated CUO. To this end, we employ multi-year canopy-top observations of micro-meteorology, ozone mixing ratios, surface energy balance components, and fluxes of ozone to derive the stomatal and non-stomatal components of the total ozone flux, combined with observation-driven ozone dry deposition simulations, using the two aforementioned representations.

2 Data and methods

2.1 Site description

Our study makes use of half-hourly observations of micro-meteorology (net radiation, air pressure, air temperature, relative humidity, precipitation, wind speed, and friction velocity) surface energy balance components and fluxes of ${\mathrm{CO}}_{\mathrm{2}}$ and ozone from two forested flux observation sites (Ispra and Hyytiälä), which are detailed below.

The Ispra forest flux station is situated in a deciduous forest in northern Italy (45.81${}^{\circ }$ N, 8.63${}^{\circ }$ E) at the European Commission Joint Research Centre (EC-JRC) in a 10 ha almost natural ecosystem mainly consisting of Quercus robur (80 %), Alnus glutinosa (10 %), Populus alba (5 %), and Carpinus betulus (3 %). Leaf area index (LAI) shows an average value of 4.1 m${}^{\mathrm{2}}$ m${}^{-\mathrm{2}}$ during the growing season . In our analysis, we rely on continuous LAI measurements unavailable at this site, which we, therefore, take from a remote sensing product derived from MODIS . The LAI range at Ispra in this product is 0.7–3.7 m${}^{\mathrm{2}}$ m${}^{-\mathrm{2}}$, scaled up to a locally measured LAI maximum of 4.5 m${}^{\mathrm{2}}$ m${}^{-\mathrm{2}}$ in July 2015 , using a seasonally varying sinusoidal scaling function. The turbulent flux measurements of surface energy balance components and ozone were performed in 2013–2015 at 36 m a.g.l. (above ground level), approximately 10 m above the canopy height of 26 m. More information regarding the measurement setup of this site can be found in .

The Hyytiälä SMEAR II (Station for Measuring Forest Ecosystem Atmosphere Relations) measurement station is located in a needleleaf forest in southern Finland (61.85${}^{\circ }$ N, 24.28${}^{\circ }$ E), with a forest cover dominated by pine trees. LAI was periodically measured at this site and varies between 2.3 and 4 m${}^{\mathrm{2}}$ m${}^{-\mathrm{2}}$. Ozone flux measurements are available for 2002–2012, with a 1-year data gap in 2006. Turbulent flux measurements are performed at 23 m a.g.l., which is 5–9 m above the forest top of 14–18 m. Ozone mixing ratios at this altitude are derived by linearly interpolating between observations at 16.8 and 33 m. More information about the measurement setup of this site and eddy covariance flux calculation can be found in and .

2.2 Observational approach

Our observational analysis, schematically depicted in Fig. a, aims to derive bulk canopy stomatal and non-stomatal resistances from canopy-top eddy covariance observations in order to estimate the magnitude of stomatal and non-stomatal ozone removal. We first derive the ozone canopy conductance ($g_{\mathrm{c},{\mathrm{O}}_{\mathrm{3}}}$) from the observed ozone dry deposition velocity ($V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$), measurement-inferred aerodynamic resistance $r_{\text{a}}$, and bulk canopy quasi-laminar layer resistance $r_{\text{b}}$ (see the Supplement).

We use the inverted Penman–Monteith equation to derive bulk canopy stomatal conductance ($g_{\text{s}}$) from canopy-top eddy covariance observations of the latent heat flux complemented with other observed variables as follows :

1$$g_{\text{s}}={\displaystyle \frac{\mathit{\lambda }E{g}_{\text{a}}\mathit{\gamma }}{\mathrm{\Delta }(R_{\text{n}}-G)+\mathit{\rho }{c}_{\text{p}}{g}_{\text{a}}\text{VPD}-\mathit{\lambda }E(\mathrm{\Delta }+\mathit{\gamma })}},$$ where $g_{\text{a}}$ is the aerodynamic conductance to water vapor (see Sect. S1 in the Supplement), $\mathit{\lambda }E$ is the latent heat flux, $\mathit{\gamma }$ is the psychrometric constant, which relates the water vapor partial pressure to air temperature, $\mathrm{\Delta }$ is the slope of the saturation vapor pressure curve, $R_{\text{n}}$ is net radiation, $G$ is the ground heat flux, $\mathit{\rho }$ is the air density, $c_{\text{p}}$ is the specific heat of air, and VPD is the vapor pressure deficit. Note that all components of Eq. () are observed or derived from observations. $g_{\text{s}}$ refers to stomatal conductance to ${\mathrm{H}}_{\mathrm{2}}\mathrm{O}$. When we refer to the stomatal conductance for ozone, we scale $g_{\text{s}}$ for the diffusivity ($D$) ratio of ozone and water vapor as follows: $g_{\mathrm{s},{\mathrm{O}}_{\mathrm{3}}}=\frac{D_{{\mathrm{O}}_{\mathrm{3}}}}{D_{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}}{g}_{\mathrm{s},{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}=\mathrm{0.61}{g}_{\mathrm{s},{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}$. Non-stomatal conductance is derived as the residual of the bulk canopy conductance and the canopy stomatal conductance, assuming that stomatal and (bulk) non-stomatal uptake are two parallel pathways (see Fig. a).

Figure 1

Schematic display representing the biophysical controls on surface ozone removal in plant canopies in the three different approaches in this study, together with their input and output variables. The combination of uptake resistances (shown as black rectangles) inside the dashed gray rectangle yields the bulk canopy resistance ($r_{\text{c}}$). In- and output variables of the mechanisms are shown in blue and green, respectively. Orange rectangles in panel (c) display the derivation of photosynthesis parameters required in MLC-CHEM; this procedure is described in more detail in Appendix . The shown resistances are the stomatal resistance ($r_{\text{s}}$), bulk canopy non-stomatal resistance ($r_{\text{ns}}$), resistance to cuticular uptake ($r_{\text{cut}}$), the resistance to in-canopy transport ($r_{\text{a,inc}}$), resistance to soil uptake ($r_{\text{soil}}$), and resistance to in-canopy transport in the upper canopy layer ($r_{\text{a,uc}}$), lower canopy layer ($r_{\text{a,lc}}$), and to the soil ($r_{\text{a,soil}}$).

[Figure omitted. See PDF]

2.3 Ozone uptake parameterizations

2.3.1 The big leaf approach

The parameterization of gaseous dry deposition in many atmospheric chemistry models is based on the resistance in the series framework introduced by , hereafter referred to as W89. The discussion below considers the implementation of the big leaf dry deposition approach in the coupled meteorology–chemistry model WRF-Chem . Other big leaf parameterizations are available with improved treatments of the stomatal e.g., and non-stomatal uptake e.g.,. However, the common use of Wesely's () parameterization in state-of-science 3D atmospheric chemistry and transport models see, e.g., motivates the choice for this scheme in our experimental setup.

Figure b depicts the resistance framework. Note that this dry deposition representation is zero-dimensional, i.e., no explicit in-canopy ozone mixing ratios are calculated. The aerodynamic resistance ($r_{\text{a}}$) is calculated following Monin–Obukhov similarity theory, and the quasi-laminar layer resistance ($r_{\text{b}}$) is estimated following . Stomatal resistance is calculated as follows :

2$$r_{\text{s}}=r_{\text{i}}\left(\mathrm{1}+{\left({\displaystyle \frac{\mathrm{200}}{R_{\text{n}}+\mathrm{0.1}}}\right)}^{\mathrm{2}}\right)\left({\displaystyle \frac{\mathrm{400}}{T_{\text{s}}(\mathrm{40}-T_{\text{s}})}}\right),$$ where $r_{\text{i}}$ (internal resistance) is a season- and land-use-dependent scaling factor, $R_{\text{n}}$ is net radiation, and $T_{\text{s}}$ is the surface temperature. $r_{\text{s}}$ is corrected for the diffusivity difference between ${\mathrm{H}}_{\mathrm{2}}\mathrm{O}$ and ozone, as explained in Sect. . In this formulation, the resistance to stomatal uptake is lowest during high-radiation conditions and for an optimum temperature of 20 ${}^{\circ }$C, reflecting that stomatal aperture follows a diurnal cycle with a peak around midday. Note that this parameterization does not explicitly account for stomatal closure due to a vapor pressure deficit or soil moisture stress. We use the non-stomatal resistances, following , which are all constant, except for the resistance to transport to the lower canopy that depends inversely on net radiation. For the soil uptake resistance, we use site-inferred values of 300 s m${}^{-\mathrm{1}}$ for Ispra and 400 s m${}^{-\mathrm{1}}$ for Hyytiälä .

2.3.2 The Multi-Layer Canopy–CHemistry Exchange Model (MLC-CHEM)

We also apply the Multi-Layer Canopy–CHemistry Exchange Model (MLC-CHEM) to evaluate simulated long-term canopy-scale ozone deposition at the two sites. This one-dimensional model explicitly simulates the canopy exchange and vertical profiles of ozone concentrations as a function of radiation, turbulent mixing, chemistry (using the Carbon Bond Mechanism, version 4 – CBM-IV), biogenic emissions following the Model for Emissions of Gases and Aerosols from Nature (MEGAN);, soil NO emissions , and (non-)stomatal uptake and their vertical gradients in the canopy. MLC-CHEM has been applied coupled to single-column and global chemistry–climate modeling studies , as well as in an offline setup for the interpretation of site-scale measurements e.g.,.

In our setup, the model consists of three layers representing the understory and the crown layer, as well as one layer aloft representing a bulk surface layer. In-canopy exchange is represented by two canopy layers whose depth depends on the canopy height ($h_{\text{c}}$), each with a layer thickness of 0.5 $h_{\text{c}}$. This two-canopy layer setup allows the simulation of in-canopy concentration and flux profiles using a computationally efficient analytical solution, allowing for the coupling of MLC-CHEM to single-column and global chemistry–climate modeling studies . Given the large gradients in radiation in the canopy, vertical profiles of radiation and radiation-dependent processes (photolysis and biogenic emissions) are calculated considering four canopy layers. The four-layer radiation profiles and biogenic emission rates are subsequently averaged over the two canopy layers for the exchange simulation. The model simulation time step is 30 min, but for processes requiring a higher temporal resolution, a sub-time-step temporal resolution is applied, which depends on the removal rate .

Micro-meteorological variables are provided as input to the model, and ozone concentrations in the upper layer are nudged to the observed above-canopy ozone concentrations to represent entrainment and advection. We use a weighting factor of 0.5, which implies that we force simulated above-canopy ozone mixing ratios to observed mixing ratios with a timescale of $\pm \mathrm{2}$ h, based on the applied temporal resolution of 0.5 h. The specific procedure to incorporate observations in our model setup is described in Sect. .

In-canopy aerodynamic resistance ($r_{\text{a}}$) is calculated as a function of canopy height, LAI, and $u_{*}$. Leaf-level stomatal conductance is calculated using the assimilation–stomatal conductance model, $A$-$g_{\text{s}}$, as follows :

3$$g_{\text{s,c,leaf}}=g_{\text{min,c}}+{\displaystyle \frac{a_{\mathrm{1}}{A}_{\text{g}}}{\left(\right[{\mathrm{CO}}_{\mathrm{2}}]-\mathrm{\Gamma })\left(\mathrm{1}+\frac{D_{\text{s}}(a_{\mathrm{1}}-\mathrm{1})}{D_{\mathrm{0}}}\right)}},$$ where $g_{\text{min,c}}$ (cuticular conductance), the constant $a_{\mathrm{1}}$, and $\mathrm{\Gamma }$ (the ${\mathrm{CO}}_{\mathrm{2}}$ compensation point) depend on the vegetation type. $A_{\text{g}}$ is gross assimilation, calculated as a function of photosynthetically active radiation (PAR), skin temperature, the internal ${\mathrm{CO}}_{\mathrm{2}}$ concentration, and the soil water content (SWC). We refer the reader to Appendix A in for more details on the calculation of $A_{\text{g}}$. $D_{\text{s}}$ is the vapor pressure deficit (VPD) leaf level, and $D_{\mathrm{0}}$ is the VPD at which stomata close. $g_{\text{s,c,leaf}}$ is calculated at the leaf level and subsequently integrated to the specific layer as a function of layer-specific LAI and PAR . This stomatal conductance representation accounts for observed increases in $g_{\text{s}}$ for an increase in ${\mathrm{CO}}_{\mathrm{2}}$ assimilation (which responds to radiation), whereas $g_{\text{s}}$ decreases as the external ${\mathrm{CO}}_{\mathrm{2}}$ concentration increases (a lower ${\mathrm{CO}}_{\mathrm{2}}$ uptake rate is needed to maintain the supply of ${\mathrm{CO}}_{\mathrm{2}}$ to the photosynthesis mechanism). $g_{\text{s}}$ also decreases as the vapor pressure deficit increases in order to minimize plant water loss through transpiration. This is a more mechanistic description of stomatal conductance compared to the big leaf approach (Eq. ), where $g_{\text{s}}$ is parameterized as a function of radiation and temperature.

The $A$-$g_{\text{s}}$ model has several degrees of freedom in determining the parameter settings. In order to derive physically appropriate settings, we tested the sensitivity of the MLC-CHEM-simulated canopy stomatal conductance ($g_{\text{s}}$) and the canopy ${\mathrm{CO}}_{\mathrm{2}}$ flux to $A$-$g_{\text{s}}$ parameter settings by comparing with observation-inferred $g_{\text{s}}$ (using Eq. ) and canopy-top $F_{{\mathrm{CO}}_{\mathrm{2}}}$ observations (see Fig. c). This procedure is described in Appendix , and the final, optimized $A$-$g_{\text{s}}$ parameters are shown in Table . With this approach, we effectively implement a realistic, observation-constrained representation canopy-top ${\mathrm{CO}}_{\mathrm{2}}$ flux and $g_{\text{s}}$ in MLC-CHEM.

Non-stomatal removal in MLC-CHEM is represented using uptake resistances taken from , , and . Analogous to W89, we adapt MLC-CHEM's default soil uptake resistance to site-inferred values of 300 s m${}^{-\mathrm{1}}$ for Ispra and 400 s m${}^{-\mathrm{1}}$ for Hyytiälä . Experimental evidence suggests increased deposition to dew-wetted leaves . MLC-CHEM accounts for this by using two distinct uptake resistances for the deposition to leaf cuticles and uptake by water films on leaves of 10${}^{\mathrm{5}}$ and 2000 s m${}^{-\mathrm{1}}$, respectively . Canopy wetness is represented by inferring the fraction of wet vegetation ($f_{\text{wet}}$) as a function of relative humidity (RH) as follows : 4$$f_{\text{wet}}=\left\{\begin{array}{ll}\mathrm{1}& \text{RH}\ge \mathrm{0.9}\\ \frac{\text{RH}-\mathrm{0.55}}{\mathrm{0.35}}& \mathrm{0.55}\le \text{RH}<\mathrm{0.9}\\ \mathrm{0}& \text{RH}<\mathrm{0.55}.\end{array}\right.$$

2.4 Experimental setup

We apply the W89 big leaf parameterization and the multi-layer ozone atmosphere–biosphere exchange parameterization to simulate total canopy ozone removal and its partitioning into stomatal and non-stomatal removal, at two locations with contrasting climate and pollution regimes, for a total of 12 site years. These simulations are compared against observation-inferred $g_{\text{s}}$ and $g_{\text{ns}}$. We restrict this analysis to daytime values (08:00–20:00 LT) during April–September, which approximately coincides with the growing season. The observational approach is known to be biased under high canopy wetness conditions due to dew formation or precipitation, and various approaches to correct for this have been reported in the literature e.g.,. We, therefore, only include data with RH $<$ 90 % and when the accumulated precipitation in the preceding 12 h is less than 0.1 mm. This set of assumptions compromises between data quality and retention of data points.

Figure 2

Time series of April–September monthly average daytime (08:00–20:00 LT) ozone dry deposition velocity for Ispra, for W89 (blue), MLC-CHEM (green), and observations (black). Solid lines and points show monthly daytime medians for simulations and observations, respectively, and shaded areas and whiskers display the interquartile range.

[Figure omitted. See PDF]

3 Results

3.1 Temporal variability in ozone dry deposition velocity

3.1.1 Monthly and interannual variability

The observed ozone uptake at Ispra is generally highest in June–August, with little interannual variability (Fig. . W89 underestimates the observed dry deposition velocity ($V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$) by $\pm \mathrm{0.1}$ cm s${}^{-\mathrm{1}}$, while MLC-CHEM reproduces the observed magnitude of $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ within 7 % in May–September. On the basis of the statistical model performance metrics in Table , there is no parameterization that consistently outperforms the other on monthly timescales. MLC-CHEM systematically overestimates ozone deposition in April. To evaluate this bias further, we performed MLC-CHEM simulations with a deactivated sink to wet leaves, motivated by the considerable uncertainty in this ozone removal pathway . This simulation resulted in the strongest decrease in $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ in the relatively humid months of April (Fig. S2), ranging from 0.15 cm s${}^{-\mathrm{1}}$ in April 2013 to 0.05 cm s${}^{-\mathrm{1}}$ in April 2015. This modification results in an improved representation of seasonality in $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$, suggesting seasonal variation in the ozone sink to wet leaves that might not be properly captured by the RH-dependent parameterization of wet leaf uptake (Eq. ).

The observed $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ at Hyytiälä is generally lower compared to Ispra, reflecting a lower leaf area and, thus, less stomatal uptake at the Finnish site. W89 and MLC-CHEM both capture the observed magnitude of $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ to within the interquartile range of observations ($\pm \mathrm{0.2}$ cm s${}^{-\mathrm{1}}$) in most years, although $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ in W89 peaks 1 month early compared to the observations. MLC-CHEM reproduces the seasonal cycle in $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ with a Pearson (temporal) correlation coefficient between simulations and observations, which is markedly higher compared to the W89 approach ($r^{\mathrm{2}}=\mathrm{0.59}$ for MLC-CHEM; $r^{\mathrm{2}}=\mathrm{0.11}$ for W89; Table ). These results suggest that MLC-CHEM better reproduces stomatal and non-stomatal removal processes, and we will investigate this further below.

The interannual variability in the ozone dry deposition velocity for Hyytiälä is 0.17 cm s${}^{-\mathrm{1}}$ and is slightly underestimated in both simulations (0.10–0.11 cm s${}^{-\mathrm{1}}$; not shown). We therefore calculated the contributions from stomatal conductance and non-stomatal conductance to the overall deposition velocity, as described in Sect. . Interannual variability in stomatal conductance is overestimated slightly by W89 and MLC-CHEM, compared to the observation-derived $g_{\text{s}}$ estimates, by 0.02 and 0.05 cm s${}^{-\mathrm{1}}$, respectively. Interannual variability in non-stomatal conductance is strongly underestimated in both simulations (0.04–0.07 cm s${}^{-\mathrm{1}}$) compared to the observed interannual variability in non-stomatal conductance (0.19 cm s${}^{-\mathrm{1}}$). The missing interannual variability in the non-stomatal deposition pathway may be due the chemical, wet leaf, and soil uptake pathways.

Figure 3

As in Fig.  but for Hyytiälä.

[Figure omitted. See PDF]

Table 1

Performance statistics for the monthly averaged simulations of $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ with W89 and MLC-CHEM (MLC). The unit is centimeters per second (cm s${}^{-\mathrm{1}}$)for mean bias error (MBE), root mean square error (RMSE), and the intercept and unitless for the other metrics. Shown are several conventionally applied performance metrics (MBE, RMSE, slope (s) and intercept (i) of a linear regression fit of simulations against observations and $r^{\mathrm{2}}$ from the ordinary least squares regression) and the index of agreement (d).

	MBE	RMSE	$r^{\mathrm{2}}$	Slope (s),	Agreement
				intercept (i)	(d)
	Ispra ($n=\mathrm{18}$ months)
W89	$-$0.09	0.14	0.24	0.79, 0.19	0.60
MLC	0.08	0.13	0.28	1.18, v0.20	0.58
	Hyytiälä ($n=\mathrm{45}$ months)
W89	$-$0.01	0.15	0.11	0.32, 0.28	0.62
MLC	0.04	0.09	0.59	1.26, $-$0.16	0.78

3.1.2 Diurnal cycles

The observed diurnal cycle of $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ at Ispra (Fig. a) is characterized by an asymmetrical pattern, with a steep morning increase that plateaus around 0.8 cm s${}^{-\mathrm{1}}$ and a decrease in the afternoon that reflects stomatal closure and reduced non-stomatal uptake. W89 underestimates the observed median daytime $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ values by $\pm \mathrm{0.1}$ cm s${}^{-\mathrm{1}}$ (20 %), while MLC-CHEM reproduces the observations within 10 %. The onset of the W89-simulated daytime $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ peak shows a 1 h time lag, with an underestimation of around $-$0.3 cm s${}^{-\mathrm{1}}$ (52 %) in the morning (06:00–10:00 LT) and an overestimation of 0.1 cm s${}^{-\mathrm{1}}$ (13 %) in the afternoon (12:00–16:00 LT). The contribution of the stomatal and non-stomatal removal in this model–observation mismatch will be discussed in Sect.. MLC-CHEM reproduces the diurnal course of $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ within 0.1 cm s${}^{-\mathrm{1}}$ throughout the day.

The observed $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ diurnal cycle at Hyytiälä (Fig. b) increases earlier during the day, compared to Ispra, and decreases later, due to the extended day length during the growing season at the Finnish site. $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ peaks at 0.5 cm s${}^{-\mathrm{1}}$ between 09:00–12:00 LT and decreases in the early afternoon due to decreasing (non-)stomatal sink ozone removal. W89 overestimates the magnitude of $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ by 0.1 cm s${}^{-\mathrm{1}}$ (22 %) in the afternoon (12:00–16:00 LT) and underestimates ozone uptake in the morning (03:00–10:00 LT) and evening (after 19:00 LT). Apart from a morning overestimation by up to 0.1 cm s${}^{-\mathrm{1}}$, MLC-CHEM reproduces the diurnal evolution of $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ well, which is apparently due to a more realistic representation of stomatal and non-stomatal removal processes.

Figure 4

Diurnal cycles of April–September ozone dry deposition velocity at Ispra (a) and Hyytiälä (b) derived from observations (black) and simulations with the W89 parameterization (blue) and MLC-CHEM (green). Lines and points show median values, and shaded areas and whiskers display the interquartile range.

[Figure omitted. See PDF]

3.2 Diurnal variability in stomatal and non-stomatal uptake

3.2.1 Ispra

Next, we analyze the stomatal and non-stomatal components of ozone deposition to further understand the model–observation agreement on diurnal timescales. Figure  shows growing season median diurnal cycles of bulk canopy conductance ($g_{\text{c}}$), canopy stomatal conductance ($g_{\text{s}}$), and non-stomatal conductance ($g_{\text{ns}}$) for Ispra in W89 and MLC-CHEM simulation and observational estimates. At Ispra, the observation-derived daytime median ozone canopy conductance is 0.87 cm s${}^{-\mathrm{1}}$ (Fig. a). The inferred daytime median stomatal conductance is as small as 0.26 cm s${}^{-\mathrm{1}}$ (gray points in Fig. b), corresponding to a daytime stomatal uptake fraction of 35 % (Fig. d). However, we found a substantial gap (of 56 %) in the energy balance closure ($Q_{\text{gap}}$), defined as the difference between net incoming radiation ($R_{\text{n}}$) and the surface energy balance components . This indicates underestimations in observed sensible and latent energy fluxes ($H$ and $LE$), which affects our observation-derived stomatal conductance. The energy balance closure issues remain after filtering the observations based on quality flags and $u_{*}$ thresholds (Fig. S1).

To resolve these energy balance closure issues, we applied a correction method that partitions $Q_{\text{gap}}$ to $H$ and $LE$ via the evaporative fraction EF $=$ $LE/(H+LE)$;. This correction increases $LE$ and $H$ by 156 and 25 W m${}^{-\mathrm{2}}$, respectively, corresponding to an evaporative fraction of 0.86. With these corrected surface energy balance components, we derive a substantially larger daytime median stomatal conductance to ozone of 0.49 cm s${}^{-\mathrm{1}}$ (black points in Fig. b, d), which is an increase of nearly 90 % with respect to the original observation-derived estimate. The $Q_{\text{gap}}$ correction also leads to a better model–observation agreement for $g_{\text{s}}$. Ozone fluxes are also affected by the surface energy balance closure gap; additional data filtering based on $u_{*}$ thresholds leads to increases in observed ozone fluxes, which exceeds 50 % in the morning and evening when absolute fluxes are low, but the effect is smaller ($<$ 15 %) during midday.

The observed diurnal cycle in canopy conductance at Ispra is better captured by MLC-CHEM compared to W89 (Fig. a). MLC-CHEM also better captures decreases in $g_{\text{c}}$ observed in the afternoon. MLC-CHEM and W89 simulate a daytime median ozone stomatal conductance of 0.43 and 0.51 cm s${}^{-\mathrm{1}}$, respectively, and, thus, agree better with the $Q_{\text{gap}}$-corrected stomatal conductance estimate derived from observations (Fig. b). The observation-derived stomatal fraction during 08:00–20:00 LT (0.62) is overestimated by W89 (0.72) and underestimated by MLC-CHEM (0.52). The observed stomatal uptake fraction increases throughout the day, from $\pm \mathrm{0.4}$ at 08:00 LT to $\pm \mathrm{0.8}$ at 18:00 LT, and this diurnal course is better reproduced by MLC-CHEM than by W89.

Observation-derived non-stomatal conductance peaks in the morning, levels off at $\pm \mathrm{0.8}$ cm s${}^{-\mathrm{1}}$ (Fig. c; gray points), and decreases in the afternoon before reaching a nighttime value of 0.1 cm s${}^{-\mathrm{1}}$. The stomatal conductance increase following $Q_{\text{gap}}$ correction leads to a reduction in the daytime average inferred non-stomatal conductance from 0.57 to 0.35 cm s${}^{-\mathrm{1}}$. This correction does, however, not affect the shape of the diurnal cycle in $g_{\text{ns}}$, which is characterized by a sharp increase in the morning and a more gradual reduction in the afternoon. Daytime non-stomatal conductance is strongly underestimated by W89 and shows little diurnal variability since most in-canopy resistances are constant and, apparently, too high. MLC-CHEM reproduces the observed diurnal evolution in non-stomatal conductance more accurately than W89 (Fig. c), apparently due to its representation of diurnal variability in processes involved in non-stomatal removal, wet leaf uptake, and in-canopy turbulence. The contributions of different removal processes to total non-stomatal uptake will be discussed in Sect. .

Figure 5

April–September median diurnal cycles of ozone bulk canopy conductance (a), canopy stomatal conductance (b), bulk non-stomatal conductance, and (c) the stomatal fraction of total ozone removal ($g_{\text{s}}{g}_{\text{c}}^{-\mathrm{1}}$; d) for Ispra. Observed medians and interquartile ranges after $Q_{\text{gap}}$ correction (OBS${}^{*}$; see text) are shown as black points and whiskers (the values prior to $Q_{\text{gap}}$ correction, denoted as OBS, are shown in gray). The median and interquartile range of W89 and MLC-CHEM are shown in blue and green, respectively. The shaded area in panel (d) highlights the nighttime period (defined as 08:00–20:00 LT) during which the stomatal flux is calculated.

[Figure omitted. See PDF]

3.2.2 Hyytiälä

At Hyytiälä, the observation-derived daytime median $g_{\text{c}}$ is 0.53 cm s${}^{-\mathrm{1}}$ (Fig. a), which is lower compared to Ispra due to lower non-stomatal ozone removal. W89 overestimates canopy conductance by up to 0.2 cm s${}^{-\mathrm{1}}$ in the afternoon, while morning and evening $g_{\text{c}}$ are underestimated. Similar to Ispra, MLC-CHEM captures the diurnal evolution in $g_{\text{c}}$ better than W89, with a peak around 09:00 LT, as in the observations, but overestimates morning canopy conductance by 0.1 cm s${}^{-\mathrm{1}}$. We did not correct for surface energy balance closure gaps for the Hyytiälä observations, since this gap was considerably smaller ($\pm \mathrm{20}$ % of $r_{\text{n}}$, without a distinct diurnal cycle) and in closer agreement to the literature-reported values for tall vegetation .

Observed stomatal conductance peaks at $\pm \mathrm{0.5}$ cm s${}^{-\mathrm{1}}$ at 10:00 LT, followed by a decrease in the afternoon (Fig. b). W89 underestimates $g_{\text{s}}$ in the morning (05:00–10:00 LT), and overestimates afternoon values by 20 %–25 %. MLC-CHEM overestimates morning stomatal conductance, but follows the observed diurnal cycle well throughout the rest of the day. The observed stomatal ozone uptake fraction is relatively constant at 0.8 (Fig. ), comparable to the upper range of stomatal uptake fraction estimates by . The stomatal fraction is well reproduced by both parameterizations, although for W89 this seems a coincidence given the misrepresented diurnal cycle in $g_{\text{c}}$ and $g_{\text{s}}$.

Observation-derived non-stomatal conductance at Hyytiälä (Fig. c) is relatively constant at 0.1 cm s${}^{-\mathrm{1}}$, except for a morning peak of 0.2 cm s${}^{-\mathrm{1}}$ around 08:00 LT that likely reflects wet leaf ozone uptake . W89 reproduces the observed daytime magnitude of $g_{\text{ns}}$ but cannot reproduce its morning peak. MLC-CHEM overestimates the nighttime non-stomatal ozone sink, in line with a study by based on a 1-month time series of ozone flux observations (August 2010), indicating that observed nighttime ozone deposition appears to reflect smaller nocturnal soil uptake efficiency than assumed. Except for an overestimation in the morning, MLC-CHEM captures the observation-inferred magnitude of non-stomatal ozone deposition well during daytime.

Figure 6

As in Fig.  but for Hyytiälä.

[Figure omitted. See PDF]

3.3 Dependence of non-stomatal deposition on driving variables

Non-stomatal ozone uptake, and its dependence on micro-meteorological and other environmental drivers, is incompletely understood. Previous studies employed statistical or process-oriented modeling to determine the contribution of driving variables to this ozone sink. In this section, we study observed and simulated relationships between the non-stomatal ozone removal fraction ($g_{\text{ns}}{g}_{\text{c}}^{-\mathrm{1}}$) and two variables (air temperature, $T_{\text{a}}$, and VPD) that we hypothesize to contribute to temporal variability in non-stomatal ozone removal. This section focuses on non-stomatal ozone removal at Ispra, since previously characterized the non-stomatal ozone sink for Hyytiälä, and we compare our findings for Ispra to their results at the end of this section.

We first determine how W89 and MLC-CHEM can reproduce the observed relationship between non-stomatal ozone removal and $T_{\text{a}}$ and VPD. We focus on the average daytime response (08:00–20:00 LT) and, subsequently, on three different periods in the diurnal cycle (06:00–20:00, 10:00–14:00, and 14:00–18:00 LT). In this manner, we can disentangle processes affecting (non-)stomatal uptake that act during different periods of the diurnal cycle (wet leaf uptake in the morning, optimal stomatal functioning during midday, and suppressed stomatal conductance during the afternoon).

The temperature response of the relative contribution of non-stomatal removal to total ozone deposition (expressed by the non-stomatal fraction, $g_{\text{ns}}{g}_{\text{c}}^{-\mathrm{1}}$) during different periods of the diurnal cycle is shown in Fig. a–d. Non-stomatal uptake decreases with temperature during the day (Fig. a). This decrease is largely driven by the morning temperature sensitivity of $g_{\text{ns}}{g}_{\text{c}}^{-\mathrm{1}}$, which shows less sensitivity to temperature later during the day (Fig. b–d). W89 underestimates the observed temperature dependence of the non-stomatal fraction throughout the day by $\pm \mathrm{0.2}$, although the morning non-stomatal fraction is higher for the lowest temperature bin (10–15 ${}^{\circ }$C). MLC-CHEM reproduces the daytime response well, as it is characterized by elevated morning non-stomatal uptake under low-temperature conditions. For most temperature bins, W89 strongly underestimates the observed variability in the non-stomatal fraction. The observed variability is also underestimated by MLC-CHEM, although to a smaller extent, and apparently indicates still missing or misrepresented deposition processes.

The observation-derived non-stomatal fraction increases with VPD during daytime (Fig. e–h), indicating that non-stomatal ozone removal decreases under dry conditions. This result contradicts an anticipated increase in the contribution by non-stomatal removal to overall canopy removal due to a VPD-induced decrease in stomatal uptake. However, the observed non-stomatal uptake also decreases in the afternoon (Fig. c) and, therefore, does not compensate for the decreasing stomatal sink with VPD. The non-stomatal fraction displays the strongest VPD sensitivity in the morning, which mainly reflects simulated wet leaf uptake under humid (i.e., low VPD) conditions. Non-stomatal removal in W89 is insensitive to VPD, and this parameterization particularly underestimates the non-stomatal fraction under humid conditions (Fig. a–b). MLC-CHEM reproduces the daytime slope between VPD and the non-stomatal fraction well.

We then perform a number of sensitivity experiments with MLC-CHEM with deactivated non-stomatal sinks to identify the role of each sink in explaining temporal variability in non-stomatal ozone removal and its dependence on $T_{\text{a}}$ and VPD. In these experiments, we exclude the contribution by wet leaf uptake, soil uptake, and in-canopy chemical removal, as well as an experiment with strongly enhanced turbulent exchange between the crown layer and the understory. In the Supplement, Sect. S2 and Fig. S6 display the results from the sensitivity analysis of non-stomatal removal at Ispra. We list the main outcomes of this section and the MLC-CHEM sensitivity analysis as follows:

• Soil deposition accounts for almost 40 % of non-stomatal removal under high-temperature conditions, reflecting a simulated increase in in-canopy turbulent transport with air temperature.

• Non-stomatal uptake is elevated under cold and humid conditions in the morning. This is consistent with MLC-CHEM-simulated wet leaf uptake, which accounts for over 20 % of the morning non-stomatal removal fraction.

• Enhanced turbulent transport from the crown layer to the understory reduces the non-stomatal uptake fraction in MLC-CHEM, as it leads to enhanced stomatal uptake in the understory.

• Chemical removal plays a minor role in the total canopy ozone sink at Ispra.

In their multivariate analysis of environmental drivers of non-stomatal ozone removal at Hyytiälä, derived that air temperature and VPD are significantly associated with variations in non-stomatal ozone removal, similar to our findings for Ispra. However, also found an explanatory role for monoterpene concentrations at Hyytiälä, while our results suggest a minor role of chemical removal at Ispra.

Figure 7

Non-stomatal ozone removal fraction $g_{\text{ns}}{g}_{\text{c}}^{-\mathrm{1}}$ binned by air temperature (a–d) and vapor pressure deficit (e–h) during June–September 2013–2015 at Ispra. Black dots and whiskers show $g_{\text{ns}}{g}_{\text{c}}^{-\mathrm{1}}$ from observations and simulations by W89 (blue points and whiskers) and MLC-CHEM (green points and whiskers). Dots and whiskers display the median and interquartile range per bin, respectively, and the number of observations in the bin is displayed at the bottom of the panels. Each column corresponds to a different time period in the diurnal cycle, namely all day (08:00–20:00 LT; a, e), morning (06:00–10:00 LT; b, f), midday (10:00–14:00 LT; c, g), and afternoon (14:00–18:00 LT; d, h). Dots and whiskers are only shown if the number of samples in the bin exceeds 10.

[Figure omitted. See PDF]

3.4 Cumulative uptake of ozone (CUO)

In the previous sections, we have shown that the seasonal evolution of ozone deposition in W89 and MLC-CHEM is relatively similar. However, there are consistent differences in daytime ozone stomatal and non-stomatal sinks between the deposition representations. In this section, we evaluate the implications of these differences in the representation of (non-)stomatal removal for determining the cumulative stomatal uptake of ozone (CUO) over the growing season, which is often used for ozone impact assessments . We here use the term CUO${}_{\text{st}}$ to refer to cumulative stomatal uptake and to distinguish this from cumulative non-stomatal ozone removal (CUO${}_{\text{ns}}$). In Fig. , we compare growing-season-integrated stomatal and non-stomatal ozone fluxes from W89 and MLC-CHEM to observation-derived estimates of total seasonal ozone uptake for both sites (Fig. a, c). Our observation-based derivation of stomatal conductance requires dry conditions (RH $<$ 90 %, and no precipitation in the preceding 12 h) to avoid overestimations in the observation-inferred stomatal conductance which lead to overestimations in CUO. However, the application of these data selection criteria also leads to a reduction in data points that hinders the calculation of CUO${}_{\text{st}}$ based on observations. In order to derive a first-order CUO${}_{\text{st}}$ estimate, we divide the cumulative stomatal uptake inferred from valid observations by the fraction of valid observations. This method serves mainly to perform a site-to-site comparison of inferred CUO${}_{\text{st}}$. Inferred CUO${}_{\text{st}}$ at Ispra varies between 61 and 72 mmol m${}^{-\mathrm{2}}$ (Fig. b). The inferred CUO${}_{\text{st}}$ in 2014 was lower compared to 2013 and 2015 due to comparatively low ozone mixing ratios, while stomatal conductance displayed less year-to-year variability. At Hyytiälä, inferred CUO${}_{\text{st}}$ varies between 39 and 41 mmol m${}^{-\mathrm{2}}$ (Fig. d), where the lower value in 2005 (29 mmol m${}^{-\mathrm{2}}$) is caused by missing data during June–August, when stomatal ozone uptake peaks. The higher inferred CUO${}_{\text{st}}$ values at Ispra compared to Hyytiälä reflect both higher stomatal conductance and ozone mixing ratios at the Italian site.

The differences between W89- and MLC-CHEM-simulated conductances are also manifested in the simulated growing season cumulative (stomatal) uptake (Fig. a–c). The cumulative total ozone flux for Ispra is underestimated by W89 ($-$10 %), while this parameterization overestimates cumulative stomatal uptake by 14 %–22 %. MLC-CHEM accurately reproduces observation-derived CUO${}_{\text{st}}$ (within 7 %) but overestimates the cumulative total flux by 15 % (Fig. a). Therefore, the model–observation agreement of the two parameterizations for the simulated cumulative total ozone removal largely reflects non-stomatal uptake differences, which deviates from observation-inferred values by $-$64 % and 51 %, respectively. At Hyytiälä, the observed cumulative total ozone flux is 15.1 mmol m${}^{-\mathrm{2}}$ and is overestimated by 34 % by W89, reflecting overestimated stomatal uptake (Fig. c). The observation-derived CUO${}_{\text{st}}$ is 12.6 mmol m${}^{-\mathrm{2}}$, compared to 15.8 mmol m${}^{-\mathrm{2}}$ in W89 ($+$28 %) and 12.3 mmol m${}^{-\mathrm{2}}$ in MLC-CHEM ($-$2.4 %). We conclude that the better representation of canopy stomatal conductance in MLC-CHEM compared to W89, particularly during the afternoon peak in ozone mixing ratios, may lead to a substantially reduced bias in the simulated growing-season-integrated (stomatal) ozone flux.

Figure 8

(a, c) Growing-season-integrated daytime (08:00–20:00 LT) and stomatal (CUO${}_{\text{st}}$; dark colors) and non-stomatal (CUO${}_{\text{ns}}$; light colors) ozone fluxes for the different years in the study period for Ispra (a) and Hyytiälä (c). Results from the W89 parameterization are shown in blue, MLC-CHEM in green, and observations in gray. Only the data points with valid observation-inferred stomatal conductance estimates are selected for this comparison; the fraction of valid data points per growing season that remains is shown at the bottom of panels (a, c). (b, d) Inferred cumulative stomatal ozone uptake (CUO) estimate at both sites (hatched bars) after dividing the season-integrated daytime stomatal ozone flux (dark gray bars) by the fraction of valid data points. Note the different $y$ axis ranges across the four panels.

[Figure omitted. See PDF]

4 Discussion

This study evaluates the potential added value of a multi-layer representation of vegetation canopies with respect to a commonly applied big leaf approach W89; for simulating ozone deposition and ozone impact metrics for forest canopies. We focus on short- to long-term temporal variability in $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ and its partitioning into stomatal and non-stomatal components, as well as the simulation of ozone impact metrics. We find that both parameterizations reasonably reproduce the observed seasonal cycle in $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$, which is in agreement with previous chemistry transport model evaluations e.g.,. Despite their comparable performance on seasonal timescales, the parameterizations deviate in their simulation of the diurnal cycle; the W89 parameterization particularly underestimates morning ozone removal by 52 % (Ispra) and 37 % (Hyytiälä), due to a combination of underestimated stomatal removal and a missing non-stomatal sink, which is likely wet leaf uptake. In the afternoon, W89 deviates less from observations at both sites ($-$13 % at Ispra; $+$22 % at Hyytiälä). Consequently, cumulative stomatal ozone uptake is overestimated by, on average, 18 % (Ispra) and 28 % (Hyytiälä) in W89 simulations, while cumulative total ozone removal deviates by $-$10 % (Ispra) and 20 % (Hyytiälä). Ozone mixing ratios typically peak in the afternoon and, thus, occur simultaneously with stomatal conductance misrepresentations, which may lead to simulated ozone fluxes overestimates using this mechanism. The multi-layer mechanism, constrained with latent energy and net ecosystem exchange (NEE) observations to optimally represent stomatal exchange, displays a better agreement with the observed ozone deposition velocity (within 10 %) and inferred cumulative stomatal and total uptake (within 15 % and 9 % for Ispra and Hyytiälä, respectively). Therefore, an accurate representation of diurnal variability in ozone uptake partitioned to stomatal and non-stomatal sinks is essential for reproducing cumulative (stomatal) ozone uptake at the land surface.

We applied a big leaf parameterization that is commonly used in (regional) atmospheric chemistry models, for example in WRF-Chem . Big leaf parameterizations advantageously depend on a limited number of routinely available meteorological variables and a simplified description of land use characteristics and can be readily applied at any location without location-specific parameter derivations . However, the empirical nature of these schemes leads to an oversimplification of in-canopy physical and chemical processes that affect atmosphere–biosphere exchange of ozone, e.g., by not accounting for stomatal closure based on the vapor pressure deficit (VPD) and soil moisture or in-canopy chemical reactions. There are big leaf versions available with a more process-based description of ozone deposition processes, particularly for stomatal conductance e.g., and non-stomatal ozone removal .

To further explore the effect of model assumptions in big leaf parameterizations, we performed a comparison between W89 and another commonly used big leaf dry deposition scheme by referred to as Z03 in Appendix . This parameterization includes a separate treatment of sunlit versus shaded leaves and the explicit treatment of water stress in the stomatal conductance calculation and includes variations in non-stomatal resistances as a function of LAI and $u_{*}$. We find that both parameterizations overestimate afternoon stomatal conductance compared to observations, while Z03 better reproduces morning $g_{\text{s}}$ (Fig. ). The differences between these parameterizations are, therefore, largely driven by differences in non-stomatal ozone removal (Fig. ). The agreement with observation-inferred non-stomatal removal depends on site-specific conditions, particularly friction velocity. Our analyses highlight potential areas of improvement in process representation that can be considered in future larger-scale modeling studies to improve the simulations of ozone deposition pathways and their temporal variability. This is particularly important for season-integrated (stomatal) ozone fluxes with big leaf parameterizations.

Our results suggest that $A_{\text{net}}$-$g_{\text{s}}$ parameterizations, as applied in MLC-CHEM, simulate stomatal conductance in good agreement with observation-inferred values throughout the diurnal cycle. Such models are sensitive to parameters typically derived at leaf level that display spatiotemporal variability. Further observational constraints on these parameters, e.g., from leaf-level ecophysiological measurements, improve the representation of stomatal conductance and biosphere–atmosphere exchange , benefitting simulations of ${\mathrm{CO}}_{\mathrm{2}}$ and ozone exchange as simulated by $A_{\text{net}}$-$g_{\text{s}}$ within MLC-CHEM. Determining these parameters from canopy-top observations is an underdetermined problem in a mathematical sense, which we circumvented by deriving a realistic set of model parameters based on a comparison with canopy-top observed NEE and observation-derived stomatal conductance, while remaining as close as possible to the original parameter set in . Choosing $A_{\text{net}}$-$g_{\text{s}}$ parameters could be formalized by applying mathematical techniques such as data assimilation .

MLC-CHEM can be driven by diagnostic variables available from the chemical transport model (CTM) output (or their driving meteorological models), favoring its implementation to represent atmosphere–biosphere fluxes of reactive compounds . In such a coupled setup, MLC-CHEM would use simulated stomatal conductance from the driving model to represent an atmosphere–biosphere exchange consistent with the model's representation of (micro-)meteorology. An implementation of $A$-$g_{\text{s}}$ with ${\mathrm{CO}}_{\mathrm{2}}$ mixing ratios, calculated online or offline, can be tested if simulated stomatal conductance estimates are unavailable.

Our analysis did not include soil moisture as a predictor of stomatal conductance. Sensitivity simulations in MLC-CHEM with observation-constrained soil water content (SWC) at different depths resulted in strong reductions in simulated NEE and $g_{\text{s}}$ during summer compared to observations, which suggests that these SWC observations are not indicative of root zone soil moisture. Nonetheless, simulations of ozone deposition and mixing ratios at various spatial scales suggest a higher predictive skill when accounting for SWC . Including this stress term is especially important in the context of projected drought risk and intensity increases in future climate scenarios that may aggravate ozone smog episodes due to a decreased stomatal sink .

Our analysis of non-stomatal ozone removal as a function of micro-meteorological drivers (air temperature and VPD) for Ispra reveals that the non-stomatal sink is elevated under low-VPD (i.e., high-RH) morning conditions, likely indicating uptake at the leaf surface in water films formed by dew . This sink is reproduced by MLC-CHEM by applying a wet canopy fraction dependent on RH and a constant wet skin uptake resistance. Observations suggest that this non-stomatal ozone sink is less important at Hyytiälä, which could be due to a lower RH threshold for the development of wet canopy conditions in MLC-CHEM compared to previous work . Since wet leaf uptake may affect simulated diurnal cycles of ozone in chemistry transport models , uptake parameterizations would benefit from better observation-based constraints on this removal process, both in terms of canopy wetness and wet leaf uptake efficiency.

Our sensitivity analysis also reveals an important role of soil deposition during the afternoon due to more active in-canopy transport. We applied a constant soil resistance to ozone uptake in our simulations, despite various environmental controls that have been identified, including air temperature, soil water content, near-surface air humidity, and soil clay content . Our results suggest a minor importance of chemical ozone removal at the two considered sites. However, we did not investigate the role of ozone scavenging by reactive sesquiterpenes or soil-emitted nitric oxide . Since most (big leaf) parameterizations work with a poorly constrained resistance to transport from the canopy-top to the soil e.g.,, the importance of the chemical and soil ozone sinks for total canopy ozone removal can be best explored with better resolved in-canopy turbulent exchange in model simulations.

We have shown that stomatal and non-stomatal sinks are not accurately reproduced using the W89 big leaf parameterization compared to observations at two forested ozone flux sites, leading to structurally biased instantaneous and growing season cumulated (stomatal) ozone flux simulations. Improved methods e.g., the DO3SE mechanism, do correct for soil moisture and VPD in the stomatal conductance calculation. Overestimated stomatal ozone fluxes also likely have implications for simulated ozone mixing ratios. Many models underestimate midday ozone mixing ratios in Europe , and a misrepresentation of land surface uptake may contribute to this bias. Therefore, an overestimated ozone deposition flux may also affect the simulation of concentration-based vegetation ozone impact metrics, such as AOT40, in the opposite direction compared to flux-based metrics. An improved model representation of the ozone deposition process will provide more confidence in the application of atmospheric chemistry models for surface air quality and vegetation ozone damage assessments.

To stimulate improvement of big leaf and multi-layer parameterizations, modelers may benefit from evaluations against existing long-term dry deposition observations in various ecosystems (e.g., forests and grassland) and for contrasting environmental conditions (e.g., during dry vs. wet seasons). Such an assessment is currently underway in stage 4 of the Air Quality Model Evaluation International Initiative AQMEII4;. Additionally, evaluation against in- and above-canopy ozone flux measurements can reveal information about non-stomatal sinks in these parameterizations, such as soil deposition and in-canopy chemical removal. Last, the application of proposed parameterizations for non-stomatal ozone sinks, such as for wet leaf uptake and soil uptake , should be tested in 3D and in single-point models of ozone deposition.

5 Conclusions

We compare ozone deposition simulations to multi-year observations at two European forested flux sites, with a focus on temporal variability, contributions from stomatal and non-stomatal sinks, and metrics for the damage incurred by ozone on vegetation. The widely used big leaf parameterization W89; and the in-canopy process-resolving MLC-CHEM model both reproduce the seasonal cycle of daytime ozone deposition velocity reasonably well, but there are important differences in the skill of the two approaches to capture the diurnal changes in ozone deposition. Specifically, W89 consistently underestimates ozone deposition velocities in the morning (by 37 %–52 %), while the afternoon model observation is somewhat smaller ($-$13 %–22 %). MLC-CHEM captures the diurnal cycle much better, with relatively small biases, in the morning ($-$9 % at Ispra, $+$17 % at Hyytiälä) and has good agreement (within 10 %) in the afternoon. Accounting for stomatal closure, wet leaf removal and in-canopy turbulent transport followed by soil uptake turns out to be important for accurately simulating ozone deposition on diurnal timescales.

The structural errors in W89 are explained by a misrepresentation of the diurnal cycle in stomatal and non-stomatal conductance. Simulations with a more recent big leaf parameterization result in similar biases regarding stomatal and non-stomatal uptake. The MLC-CHEM model, constrained by local observations of diurnal ${\mathrm{CO}}_{\mathrm{2}}$ and latent energy fluxes, captures stomatal and non-stomatal ozone conductance better. As a result, W89 systematically overestimates cumulative ozone uptake by 20 %–30 % in the growing season at Ispra and Hyytiälä, whereas MLC-CHEM reproduces cumulative ozone uptake within 3 % at both sites. We conclude that MLC-CHEM, nudged with observation-inferred stomatal conductance, accurately describes non-stomatal uptake processes and vegetation ozone impact metrics.

Sensitivity tests with MLC-CHEM for Ispra point out that, in relatively cold and humid conditions, ozone deposition on wet leaves appears to explain up to 20 % of the non-stomatal ozone sink. During high-temperature conditions characterized by efficient in-canopy transport, enhanced uptake by soils accounts for up to 40 % of non-stomatal ozone deposition. The tests suggest a minor role for chemical destruction of ozone at Ispra.

Our results indicate that current model representations of stomatal and non-stomatal ozone uptake by vegetation, often based on W89, should be thoroughly evaluated. This study provides a strategy for such evaluations and shows how a more detailed, canopy-resolving model driven by ancillary measurements of ${\mathrm{CO}}_{\mathrm{2}}$ and energy fluxes can provide more realistic estimates of ozone deposition and vegetation ozone impact metrics.

Appendix A

$A$-$g_{\text{s}}$ optimization

Prior to applying MLC-CHEM to analyze ozone fluxes at our study sites, we first paid attention to simulations of the canopy ${\mathrm{CO}}_{\mathrm{2}}$ flux ($F_{{\mathrm{CO}}_{\mathrm{2}}}$) and canopy stomatal conductance ($g_{\text{s}}$) to ensure that the photosynthesis parameterization ($A$-$g_{\text{s}}$) functions satisfactorily. An initial simulation with the default settings for the C${}_{\mathrm{3}}$ vegetation class resulted in a strongly overestimated $F_{{\mathrm{CO}}_{\mathrm{2}}}$ compared to observations at both sites (see Table ). This is accompanied by strong overestimation of the canopy stomatal conductance at Ispra, while MLC-CHEM slightly underestimates stomatal conductance at Hyytiälä.

The default $A$-$g_{\text{s}}$ settings were derived for low vegetation such as grassland and crops and are therefore not necessarily representative for forest canopies. We performed a sensitivity analysis of simulated $F_{{\mathrm{CO}}_{\mathrm{2}}}$ and $g_{\text{s}}$ to $A$-$g_{\text{s}}$ model parameters in order to determine optimized parameter sets for our simulations. These settings are given in Table . We found that strongly overestimated $F_{{\mathrm{CO}}_{\mathrm{2}}}$ is largely caused by a presumed high reference mesophyll conductance ($g_{\mathrm{m},\mathrm{298}}$), leading to overestimated transport of ${\mathrm{CO}}_{\mathrm{2}}$ in the plant's chloroplast. Our reductions in $g_{\mathrm{m},\mathrm{298}}$ are in better correspondence with previously reported estimates of 0.8–2.0 mm s${}^{-\mathrm{1}}$ for different forest plant functional types . At Ispra, we additionally modified the mesophyll conductance temperature response curve, which differs between plant species , to improve the amplitude of the seasonal cycle in simulated $F_{{\mathrm{CO}}_{\mathrm{2}}}$. At Hyytiälä, the maximum internal ${\mathrm{CO}}_{\mathrm{2}}$ concentration ($f_{\mathrm{0}}$, given as a fraction of the external ${\mathrm{CO}}_{\mathrm{2}}$ concentration) was increased to improve the correspondence with observation-derived $g_{\text{s}}$.

Our observational constraints to $A$-$g_{\text{s}}$ lead to improved simulations of $g_{\text{s}}$ and $F_{{\mathrm{CO}}_{\mathrm{2}}}$ (Table ). The parameter changes additionally affect the simulation of the ozone dry deposition velocity ($V_{\text{d}}$(${\mathrm{O}}_{\mathrm{3}}$)), as shown in Table . At Ispra, the strong reduction in stomatal conductance leads to an underestimation in $V_{\text{d}}$(O${}_{\mathrm{3}}$) (MBE $=$ $-$0.12 cm s${}^{-\mathrm{1}}$), while the other statistical metrics indicate a modest model improvement. At Hyytiälä, the growing season model overestimation is slightly reduced from 0.04 to 0.02 cm s${}^{-\mathrm{1}}$. Our approach results in a reduced model bias at the two study sites, particularly for $F_{{\mathrm{CO}}_{\mathrm{2}}}$, while taking care to stay as close as possible to the original parameter set.

Table A1

$A$-$g_{\text{s}}$ parameter settings used in MLC-CHEM simulations. The first column indicates the default C${}_{\mathrm{3}}$ settings from , and the other two columns show the optimal settings from our analysis. A dash (–) indicates that a parameter is unchanged with respect to the default C${}_{\mathrm{3}}$ value.

	C${}_{\mathrm{3}}$ (reference)	Ispra	Hyytiälä
$g_{\mathrm{m},\mathrm{298}}$ (mm s${}^{-\mathrm{1}}$)	7.0	1.5	1.5
$f_{\mathrm{0}}$ (–)	0.89	–	0.99
$g_{\text{m,T1}}$ (K)	278	283	–
$g_{\text{m,T2}}$ (K)	301	306	–
$A_{\text{m,max,T1}}$ (K)	281	286	–

Table A2

Model performance statistics of MLC-CHEM before and after $A$-$g_{\text{s}}$ optimization for canopy stomatal conductance and CO${}_{\mathrm{2}}$ flux. Shown are several conventionally applied performance metrics (MBE, RMSE, slope (s), and intercept (i) of a linear regression fit of simulations against observations and $r^{\mathrm{2}}$ from ordinary least squares regression), as well as the index of agreement (d). The units are centimeters per second (cm s${}^{-\mathrm{1}}$) and micromoles per square meter per second ($\mathrm{\mu }$mol m${}^{-\mathrm{2}}$ s${}^{-\mathrm{1}}$), respectively, unless indicated otherwise.

	Hyytiälä				Ispra
	$g_{\text{s}}$ (cm s${}^{-\mathrm{1}}$)		$F_{{\mathrm{CO}}_{\mathrm{2}}}$ ($\mathrm{\mu }$mol m${}^{-\mathrm{2}}$ s${}^{-\mathrm{1}}$)		$g_{\text{s}}$ (cm s${}^{-\mathrm{1}}$)		$F_{{\mathrm{CO}}_{\mathrm{2}}}$ ($\mathrm{\mu }$mol m${}^{-\mathrm{2}}$ s${}^{-\mathrm{1}}$)
	REF	OPT	REF	OPT	REF	OPT	REF	OPT
MBE	$-$0.04	$-$0.07	$-$20.1	$-$5.1	0.56	0.26	$-$45.8	$-$17.0
RMSE	0.43	0.41	21.1	6.1	0.89	0.66	50.3	18.4
$r^{\mathrm{2}}$ (–)	0.12	0.22	0.64	0.46	0.16	0.16	0.57	0.61
s (–), i	0.77, 0.17	1.01, 0.07	0.36, 2.34	0.65, 0.62	0.42, 0.05	0.69, $-$0.02	0.25, 3.15	0.62, 6.50
d (–)	0.48	0.56	0.30	0.63	0.49	0.52	0.29	0.55

Table A3

As in Table  but for $V_{\text{d}}$(O${}_{\mathrm{3}}$) (in centimeters per second; cm s${}^{-\mathrm{1}}$).

	Hyytiälä		Ispra
	REF	MOD	REF	MOD
MBE	0.04	0.02	0.11	$-$0.01
RMSE	0.17	0.17	0.47	0.32
$r^{\mathrm{2}}$ (–)	0.37	0.39	0.28	0.45
s (–), i	1.93, $-$0.01	0.86, 0.04	0.45, 0.32	0.89, 0.08
d (–)	0.73	0.76	0.69	0.80

Appendix B Comparison between two big leaf parameterizations

In order to derive more generic conclusions about big leaf parameterizations, we considered another commonly applied parameterization that has recently been extended to different gases by . This big leaf formulation (hereafter Z03) differs compared to the parameterization (hereafter W89) in several aspects, as follows: (1) Z03 calculates stomatal conductance for sunlit and shaded leaves differently, (2) stomatal conductance is affected by VPD and soil moisture stress, and (3) non-stomatal resistances contain seasonal and diurnal variability due to dependencies on leaf area index and friction velocity ($u_{*}$). This model version was derived from with two modifications. First, we adapted the soil resistance to locally derived values of 400 s m${}^{-\mathrm{1}}$ (Hyytiälä) and 300 s m${}^{-\mathrm{1}}$ (Ispra), similar to W89 and MLC-CHEM (see Sect. 2). The implementation by relies on observed canopy wetness, which is not available for our two study sites. We therefore parameterize canopy wetness as a function of relative humidity, analogous to MLC-CHEM (Eq. ). In this section, we compare simulations by W89 and Z03 to observations of the ozone dry deposition velocity and observation-inferred stomatal and non-stomatal conductance.

Figure  shows multi-year growing season median diurnal cycles of $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$, $g_{\text{s}}$ and $g_{\text{ns}}$ for Ispra and Hyytiälä. From this analysis, we conclude that W89 and Z03 perform similarly for Ispra compared against observed $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ (Fig. a). Z03 better captures the early morning onset of $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ for Hyytiälä than W89, but more strongly overestimates midday and afternoon $V_{\text{d}}\left({\mathrm{O}}_{\mathrm{3}}\right)$ compared to observations (Fig. b). Both parameterizations overestimate midday and afternoon $g_{\text{s}}$, while Z03 better captures the observed morning and afternoon $g_{\text{s}}$ values than W89 (Fig. c, d). For $g_{\text{ns}}$, there is no parameterization that performs best for the two sites. Both parameterizations underestimate observation-inferred $g_{\text{ns}}$ at Ispra (corrected for energy balance closure gaps; see Sect. ), while W89 better captures the magnitude of observation-inferred $g_{\text{ns}}$ (although Z03 better reproduces the shape of the diurnal cycle). This suggests that the $g_{\text{ns}}$ dependence on $u_{*}$ is less strong in the observations than is suggested in the Z03 parameterization; i.e., a sensitivity experiment with doubled $u_{*}$ values for Ispra results in daytime $g_{\text{ns}}$ values of 0.2–0.35 cm s${}^{-\mathrm{1}}$, which is an increase by a factor 2.3–2.8. Based on our findings, we conclude that the different representation of non-stomatal ozone removal drives the differences between W89 and Z03, but the magnitude of these differences depends on site-specific conditions.

Figure B1

Comparison of the dry deposition parameterizations W89 and Z03 against the observed dry deposition velocity (a, b) and observation-inferred stomatal conductance (c, d) and non-stomatal conductance (e, f) for Ispra and Hyytiälä (left and right panels, respectively). Lines and shaded areas (points and whiskers) show April–September median and interquartile range of the simulations (observations).

[Figure omitted. See PDF]

Code and data availability

MLC-CHEM source code and model output are available upon request to the corresponding author.

The supplement related to this article is available online at: https://doi.org/10.5194/acp-21-18393-2021-supplement.

Author contributions

AJV, LNG, and KFB designed the experiments, and AJV performed the simulations and data analysis. IG, IM, and GM provided observational data and expertise on the use thereof in this study. AJV wrote the paper, with contributions from all co-authors.

Competing interests

At least one of the (co-)authors is a member of the editorial board of Atmospheric Chemistry and Physics. The peer-review process was guided by an independent editor, and the authors also have no other competing interests to declare.

Disclaimer

Publisher’s note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Acknowledgements

The authors acknowledge the public availability of observational data sets used in this study via the European Fluxes Database Cluster (http://gaia.agraria.unitus.it, last access: 12 July 2021) and the Finnish research data portal (https://etsin.fairdata.fi, last access: 12 July 2021).

Financial support

This research has been supported by the Nederlandse Organisatie voor Wetenschappelijk Onderzoek (grant no. ALW-GO 16/17). The observational infrastructure at Hyytiälä is funded under the ACCC flagship by the Academy of Finland (grant no. 337549).

Review statement

This paper was edited by Leiming Zhang and reviewed by two anonymous referees.