Separating above-canopy CO₂ and O₂ measurements into their atmospheric and biospheric signatures

Atmospheric tracers are often used to interpret the local CO₂ budget, where measurements at a single height are assumed to represent local flux signatures. Alternatively, these signatures can be derived from direct flux measurements or by using fluxes derived from measurements at multiple heights. In this study, we contrast interpretation of surface CO₂ exchange from tracer measurements at a single height to measurements at multiple heights. Specifically, we analyse the ratio between atmospheric O₂ and CO₂ (exchange ratio, ER) above a forest. We consider the following two alternative approaches: the exchange ratio of the forest (ER_forest) obtained from the ratio of the surface fluxes of O₂ and CO₂ derived from measurements at multiple heights, and the exchange ratio of the atmosphere (ER_atmos) obtained from changes in the O₂ and CO₂ mole fractions over time measured at a single height. We investigate the diurnal cycle of both ER signals to better understand the biophysical meaning of the ER_atmos signal. We have combined CO₂ and O₂ measurements from Hyytiälä, Finland, during spring and summer of 2018 and 2019 with a conceptual land–atmosphere model to investigate the behaviour of ER_atmos and ER_forest. We show that the CO₂ and O₂ signals as well as their resulting ERs are influenced by climate conditions such as variations in soil moisture and temperature, for example during the 2018 heatwave. We furthermore show that the ER_atmos signal obtained from single-height measurements rarely represents the forest exchange directly, mainly because it is influenced by entrainment of air from the free troposphere into the atmospheric boundary layer. The influence of these larger-scale processes can lead to very high ER_atmos values (even larger than 2), especially in the early morning. These high values do not directly represent carbon cycle processes, but are rather a mixture of different signals. We conclude that the ER_atmos signal provides only a weak constraint on local-scale surface CO₂ exchange, and that ER_forest above the canopy should be used instead. Single-height measurements always require careful selection of the time of day and should be combined with atmospheric modelling to yield a meaningful representation of forest carbon exchange. More generally, we recommend always measuring at multiple heights when using multi-tracer measurements to study surface CO₂ exchange.

1 Introduction

Rising atmospheric CO₂ levels resulting from fossil fuel combustion and land-use change emissions, which are moderated by uptake by the terrestrial biosphere and oceans, require comprehensive assessment of carbon exchange at local and global scales (Friedlingstein et al., 2022). Atmospheric O₂ serves as a valuable tracer, enhancing our understanding of carbon exchange due to the close linkage between O₂ and CO₂ in carbon cycle processes such as fossil fuel combustion, photosynthesis, and respiration (Manning and Keeling, 2006; Worrall et al., 2013; Keeling and Manning, 2014; Bloom, 2015; Hilman et al., 2022). The exchange ratio (ER $=-{\mathrm{O}}_{\mathrm{2}}/{\mathrm{CO}}_{\mathrm{2}}$), denoted as the number of moles of O₂ exchanged per mole of CO₂, represents the specific link between O₂ and CO₂ for different processes (Keeling et al., 1998). Long-term O₂ and CO₂ measurements allow us to derive the global ocean carbon sink (Stephens et al., 1998; Rödenbeck et al., 2008; Tohjima et al., 2019) and to estimate changes in fossil fuel emissions (Pickers et al., 2022; Ishidoya et al., 2020; Rödenbeck et al., 2023).

Figure 1 Schematic overview of diurnal cycles of the surface fluxes and mole fractions of atmospheric O₂ and CO₂ above a forest canopy. The figure illustrates the dominant processes throughout the day, with forest exchange dominating the nocturnal and afternoon periods, while early morning signals are primarily influenced by entrainment of air from the residual layer or the free troposphere. The surface fluxes of O₂ and CO₂ result in the exchange ratio signal of the forest (ER_forest), while the changes in the mole fractions of O₂ and CO₂ over time can lead to variations of the exchange ratio signal of the atmosphere (ER_atmos). Note that the term “surface fluxes” refers to the fluxes from the surface layer, which includes the vegetation layer including the top of the canopy. The surface layer is the lowest 10 % of the boundary layer, where the surface directly influences the atmospheric boundary layer.

For global applications, a constant ER of 1.1 (mol mol⁻¹) is assumed for the terrestrial biosphere (Severinghaus, 1995). However, the ER of terrestrial biosphere exchange is not uniform at smaller scales; it varies between ecosystems and over time (Angert et al., 2015; Bloom, 2015; Battle et al., 2019; Hilman et al., 2022). Measuring the ERs of ecosystems and the underlying gross processes facilitates the partitioning of net ecosystem exchange (NEE) into gross primary production (GPP) and total ecosystem respiration (TER) (Ishidoya et al., 2015; Faassen et al., 2023), which is still challenging (Reichstein et al., 2005). The ER for net ecosystem exchange can be determined from the ratio of the net turbulent surface fluxes of O₂ and CO₂ above the canopy, referred to as ER_forest (see Fig. 1). The O₂ surface fluxes can be inferred from the vertical gradient: the difference between O₂ mole fraction measurements at multiple heights, together with a turbulent exchange coefficient. Currently, available instruments do not allow eddy covariance (EC) O₂ measurements. The ER_forest signal predominantly represents forest exchange occurring in and below the canopy (small-scale processes), comprising the individual ERs of TER (ER_r) and GPP (ER_a) (Ishidoya et al., 2013, 2015; Faassen et al., 2023). Alternatively, the net ecosystem ER has been estimated based on measurements of O₂ and CO₂ mole fractions in the atmosphere at a single height above the canopy. This is referred to as ER_atmos (Fig. 1) and is defined as the change in O₂ and CO₂ mole fractions over time (Seibt et al., 2004; Battle et al., 2019; Faassen et al., 2023).

In our recent study (Faassen et al., 2023), we presented a comprehensive comparison of the diurnal behaviour of ER_forest and ER_atmos using measurements collected above a boreal forest in Hyytiälä, Finland. Our analysis revealed that during the afternoon (the photosynthesis-dominant period in Fig. 1), the ER_atmos signal approaches the ER_forest value, although they did not converge completely. Furthermore, we showed that during the entrainment-dominant period (see Fig. 1), the ER_atmos signal strongly exceeded the expected ER value for biosphere exchange, which is typically around 1.1 (Severinghaus, 1995), and even surpassed 2.0. Such high ER values (>2.0) cannot be attributed to a single process such as photosynthesis, respiration, or fossil fuel combustion, as their ER values are below 2.0. We proposed that the high ER_atmos signal was likely influenced by large-scale processes, specifically the entrainment of air from the free troposphere into the boundary layer (Faassen et al., 2023). Seibt et al. (2004) and Yan et al. (2023) also argue that ER_atmos cannot capture the ER signal of a forest. In contrast, in the studies by Ishidoya et al. (2013, 2015), ER_forest and ER_atmos do result in similar values when small-scale processes dominate over large-scale processes. In Faassen et al. (2023), we concluded that an atmospheric model was needed to interpret the observed diurnal signals of ER_atmos and ER_forest. The current study delivers this model-based analysis.

Until now, atmospheric O₂ above forest canopies has primarily been modelled with relatively simple one-box models that use only the surface components, lacking implementation of boundary layer dynamics such as entrainment and boundary layer growth (Seibt et al., 2004; Ishidoya et al., 2013). Understanding how mole fractions and, consequently, how ER_atmos evolves throughout the day requires accounting for these critical processes. Yan et al. (2023) recently modelled O₂ and CO₂ within and below a canopy using a multi-layer model and showed that ER_atmos and ER_forest have diurnal and annual patterns. However, ER_atmos was treated as a constant value above the canopy and boundary layer dynamics were not accounted for. To expand on the work by Yan et al. (2023) and gain further insight into the diurnal ER_atmos behaviour above a canopy, in this study, we use the mixed-layer Chemistry Land-surface Atmosphere Soil Slab (CLASS) model (Vilà-Guerau de Arellano et al., 2015). In short, the model is able to represent the thermodynamics and biophysical processes associated with the diurnal variation in the boundary layer and can provide insights into the processes contributing to ER_atmos formation. Additionally, the model facilitates the analysis of ER_atmos behaviour under more extreme conditions such as droughts or heatwaves.

In this study, we aim to enhance our understanding of single-height O₂ and CO₂ measurements and the resulting ER_atmos signal, as observed above the canopy, and we propose a new relationship between the ER_atmos and ER_forest signal. We seek to determine whether single-height O₂ and CO₂ measurements can be employed to estimate the ecosystem's ER despite the aforementioned limitations. Additionally, we explore whether the ER_atmos signal constrains boundary layer dynamics, and we identify cases where large-scale processes (e.g. entrainment of background air) influence the signal of small-scale processes (e.g. NEE) by analysing different diurnal regimes of ER_forest and ER_atmos. We combine measurements from campaigns in Hyytiälä, Finland, during spring and summer of 2018 and 2019 with an analysis of the mixed-layer CLASS model. This combined approach allows us to address the following research questions: (1) when does ER_atmos represent local forest exchange processes and become equal to ER_forest and (2) what is the underlying physical explanation for the high ER_atmos values observed in the recent study by Faassen et al. (2023)?

In this paper we first derive a theoretical relationship between ER_atmos and ER_forest that can help us to understand which components influence the diurnal cycle of ER_atmos and when ER_atmos should indicate the same processes as ER_forest (Sect. 2). To evaluate the diurnal cycle of ER_atmos we combine observational data with the CLASS model (Sect. 3). We then show the model evaluation and the ER_atmos and ER_forest model results in Sect. 4.2, and we analyse different cases to explain the diurnal behaviour of ER_atmos during distinct periods of the day and investigate when ER_atmos represents forest exchange (Sect. 4.3). Next, we place our results in perspective and show how ER_atmos should (not) be used (Sect. 5). Finally, we present our conclusions on the physical explanations for the differences between the diurnal behaviour of both ER_atmos and ER_forest.

Figure 2 Vertical profiles of potential temperature (θ) measured by radiosondes at Hyytiälä on 12 July 2019 (a) and 24 July 2018 (b). The observations are conceptualized (black lines) to show (1) the well-mixed profiles at different time steps, (2) the jumps between the boundary layer and the free troposphere, and (3) the lapse rate in the free troposphere; 1, 2, and 3 are used to initialize the CLASS model. Panel (c) gives the theoretical vertical profiles of O₂ and CO₂ for the early morning (M) and late afternoon (A). The sizes of the arrows indicate the effects of entrainment (dashed lines) and the surface fluxes (solid lines) on the vertical profiles.

2 Fundamental concepts

2.1 The mixed-layer theory

The CLASS land–atmosphere model (Vilà-Guerau de Arellano et al., 2015) is based on the mixed-layer theory which assumes that scalars (such as O₂, CO₂, θ) are constant with height in the atmospheric boundary layer (Lilly, 1968; Tennekes, 1973). Figure 2 illustrates these assumptions for potential temperature (θ), O₂, and CO₂. Within the mixed-layer theory, no distinct surface layer exists, and a capping inversion links the mixed-layer value (the bulk constant value) with the lapse rate of the free troposphere. This inversion, termed the “jump” (Δ_(ft-bl)), represents the difference of a scalar (e.g. the CO₂ mole fraction) between the atmospheric boundary layer and the free troposphere. The free troposphere is represented by a linear change of the scalar with height (the lapse rate).

CLASS describes the well-mixed layer with a scalar constant in height (Fig. 2). This scalar (ϕ) can then be solved in the mixed layer with the following equation (Vilà-Guerau de Arellano et al., 2015):

$$\begin{array}{}\text{(1)}& {\displaystyle \frac{\partial \mathit{\varphi }}{\partial t}}={\displaystyle \frac{(\stackrel{\mathrm{‾}}{w^{\prime }{\mathit{\varphi }}^{\prime }}{)}_{\mathrm{s}}-(\stackrel{\mathrm{‾}}{w^{\prime }{\mathit{\varphi }}^{\prime }}{)}_{\mathrm{e}}}{h}}-\mathrm{adv}\left(\mathit{\varphi }\right),\end{array}$$

where ∂ϕ/∂t is the tendency (i.e. change over time) of a generic well-mixed scalar, w^′ are the deviations of the mean for w which is the vertical wind speed, and ϕ^′ are the deviations from the mean for a scalar ϕ. The term $(\stackrel{\mathrm{‾}}{w^{\prime }{\mathit{\varphi }}^{\prime }}{)}_{\mathrm{s}}$ is the surface flux of ϕ and represents the small-scale processes, $(\stackrel{\mathrm{‾}}{w^{\prime }{\mathit{\varphi }}^{\prime }}{)}_{\mathrm{e}}$ is the entrainment flux, h is the boundary layer height, and adv(ϕ) is the horizontal advection of scalar ϕ into the well-mixed layer. In contrast to the local surface exchange, $(\stackrel{\mathrm{‾}}{w^{\prime }{\mathit{\varphi }}^{\prime }}{)}_{\mathrm{e}}$ and adv(ϕ) represent large-scale processes.

The entrainment flux is dependent on the entrainment velocity and the jump:

$$\begin{array}{}\text{(2)}& (\stackrel{\mathrm{‾}}{w^{\prime }{\mathit{\varphi }}^{\prime }}{)}_{\mathrm{e}}=-w_{\mathrm{e}}\cdot {\mathrm{\Delta }}_{\text{(ft-bl)}}\mathit{\varphi }=\left({\displaystyle \frac{\partial h}{\partial t}}-w_{\mathrm{sub}}\right)\cdot {\mathrm{\Delta }}_{\text{(ft-bl)}}\mathit{\varphi },\end{array}$$

where w_e is the entrainment velocity, Δ_(ft-bl)ϕ is the jump between the free troposphere and the atmospheric boundary layer, and w_sub is the mean vertical subsidence velocity normally associated with high-pressure systems. We assume w_sub to be negligible, because our focus does not lie on the influence of synoptic scale processes.

Δ_(ft-bl)ϕ changes over time (see Fig. 2) and depends on the surface fluxes and the air that is entrained from the free troposphere (see Eq. 1):

$$\begin{array}{}\text{(3)}& {\displaystyle \frac{\partial {\mathrm{\Delta }}_{\text{(ft-bl)}}\mathit{\varphi }}{\partial t}}={\mathit{\gamma }}_{\mathit{\varphi }}\cdot w_{\mathrm{e}}-{\displaystyle \frac{\partial \mathit{\varphi }}{\partial t}},\end{array}$$

where γ_ϕ is the lapse rate of ϕ in the free troposphere and $\partial \mathit{\varphi }/\partial t$ is the change over time of the well-mixed scalar ϕ (i.e. in the boundary layer).

Lastly, the growth of the boundary layer height $\left(\frac{\partial h}{\partial t}\right)$ effectively determines the entrainment velocity and therefore the entrainment flux of a certain scalar. The growth of the boundary layer is caused by the virtual potential temperature (θ_v), also called buoyancy:

$$\begin{array}{}\text{(4)}& {\displaystyle \frac{\partial h}{\partial t}}=-{\displaystyle \frac{(\stackrel{\mathrm{‾}}{w^{\prime }{\mathit{\theta }}_{\mathrm{v}}^{\prime }}{)}_{\mathrm{e}}}{{\mathrm{\Delta }}_{\text{(ft-bl)}}{\mathit{\theta }}_{\mathrm{v}}}}+w_{\mathrm{sub}},\end{array}$$

where θ_v is the virtual potential temperature (i.e. potential temperature of dry air) and w_sub is the subsidence velocity. For more details on these equations, see Vilà-Guerau de Arellano et al. (2015) and Sect. 3.2.2 and Appendix A2 for the application of O₂.

2.2 Theoretical relationship between ER_atmos and ER_forest

The ER signal of the forest (ER_forest) is defined as (Faassen et al., 2023)

$$\begin{array}{}\text{(5)}& {\mathrm{ER}}_{\mathrm{forest}}=-{\displaystyle \frac{(F_{{\mathrm{O}}_{\mathrm{2}}}{)}_{\mathrm{s}}}{(F_{{\mathrm{CO}}_{\mathrm{2}}}{)}_{\mathrm{s}}}}\approx -{\displaystyle \frac{-K_{\mathit{\varphi }}\cdot {\mathrm{\Delta }}_{\left(z\right)}{\mathrm{O}}_{\mathrm{2}}/\mathrm{\Delta }z}{-K_{\mathit{\varphi }}\cdot {\mathrm{\Delta }}_{\left(z\right)}{\mathrm{CO}}_{\mathrm{2}}/\mathrm{\Delta }z}},\end{array}$$

where $(F_{{\mathrm{O}}_{\mathrm{2}}}{)}_{\mathrm{s}}$ and $(F_{{\mathrm{CO}}_{\mathrm{2}}}{)}_{\mathrm{s}}$ are the mean net turbulent surface fluxes of O₂ and CO₂, respectively, over a certain time period above the canopy and can be derived from the vertical gradient of O₂ (Δ_(z)O₂) and CO₂ (Δ_(z)CO₂) measurements at two heights together with an exchange coefficient following the K-theory (K_ϕ) (Faassen et al., 2023). Note that the K-theory does not apply when one of the measurement levels is inside the canopy. For readability, we write the surface fluxes for both O₂ and CO₂ as F_ϕ, instead of $(\stackrel{\mathrm{‾}}{w^{\prime }{\mathit{\varphi }}^{\prime }}{)}_{\mathrm{s}}$ that was used above for the general theory.

The ER signal of the atmosphere (ER_atmos) is defined as (Faassen et al., 2023)

$$\begin{array}{}\text{(6)}& {\mathrm{ER}}_{\mathrm{atmos}}=-{\displaystyle \frac{\partial {\mathrm{O}}_{\mathrm{2}}/\partial t}{\partial {\mathrm{CO}}_{\mathrm{2}}/\partial t}}\approx -{\displaystyle \frac{{\mathrm{\Delta }}_{\left(t\right)}{\mathrm{O}}_{\mathrm{2}}}{{\mathrm{\Delta }}_{\left(t\right)}{\mathrm{CO}}_{\mathrm{2}}}},\end{array}$$

where Δ_(t)O₂ and Δ_(t)CO₂ are the changes in the O₂ and CO₂ mole fractions over time (tendencies) at a single height. Linear regression between O₂ and CO₂ can be applied, and the slope gives the ER_atmos value for a certain event or time period. For this study, linear regression was applied for the three periods described in Sect. 4.2.1 for the observations (1 value per 30 min) and the CLASS model output (1 value per 10 s).

According to the mixed-layer theory described above, the tendencies in Eq. (6) depend on the surface and entrainment fluxes, together with the boundary layer height (h) (see Eq. 1). Equation (6) can be rewritten by implementing Eq. (1):

$$\begin{array}{}\text{(7)}& {\mathrm{ER}}_{\mathrm{atmos}}=-{\displaystyle \frac{\left(\right(F_{{\mathrm{O}}_{\mathrm{2}}}{)}_{\mathrm{s}}-\left(F_{{\mathrm{O}}_{\mathrm{2}}}{)}_{\mathrm{e}}\right)/h}{\left(\right(F_{{\mathrm{CO}}_{\mathrm{2}}}{)}_{\mathrm{s}}-\left(F_{{\mathrm{CO}}_{\mathrm{2}}}{)}_{\mathrm{e}}\right)/h}},\end{array}$$

where $(F_{{\mathrm{O}}_{\mathrm{2}}}{)}_{\mathrm{s}}$ and $(F_{{\mathrm{CO}}_{\mathrm{2}}}{)}_{\mathrm{s}}$ are the net surface fluxes of O₂ and CO₂, and $(F_{{\mathrm{O}}_{\mathrm{2}}}{)}_{\mathrm{e}}$ and $(F_{{\mathrm{CO}}_{\mathrm{2}}}{)}_{\mathrm{e}}$ are the entrainment fluxes of O₂ and CO₂, respectively. For simplicity, we ignored the advection term in Eq. (1) here, but we will add it later (Eq. 9). As shown in Eq. (2), the entrainment flux depends on the entrainment velocity (w_e) and the jump between the free troposphere and the boundary layer (Δ_(ft-bl)ϕ). Combining the definition of ER_forest (Eq. 5) with Eq. (2) allows us to rewrite Eq. (7) as

$$\begin{array}{}\text{(8)}& \begin{array}{rl}{\mathrm{ER}}_{\mathrm{atmos}}& ={\mathrm{ER}}_{\mathrm{forest}}\cdot \left({\displaystyle \frac{\mathrm{1}+\frac{w_{\mathrm{e}}\cdot {\mathrm{\Delta }}_{\text{(ft-bl)}}{\mathrm{O}}_{\mathrm{2}}}{(F_{{\mathrm{O}}_{\mathrm{2}}}{)}_{\mathrm{s}}}}{\mathrm{1}+\frac{w_{\mathrm{e}}\cdot {\mathrm{\Delta }}_{\text{(ft-bl)}}{\mathrm{CO}}_{\mathrm{2}}}{(F_{{\mathrm{CO}}_{\mathrm{2}}}{)}_{\mathrm{s}}}}}\right)\\ & ={\mathrm{ER}}_{\mathrm{forest}}\cdot \left({\displaystyle \frac{\mathrm{1}+{\mathit{\beta }}_{{\mathrm{O}}_{\mathrm{2}}}}{\mathrm{1}+{\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}}}\right),\end{array}\end{array}$$

where Δ_(ft-bl)O₂ and Δ_(ft-bl)CO₂ are the jumps of O₂ and CO₂ between the free troposphere and the boundary layer, and β_ϕ is the ratio between the entrainment flux and the surface flux (Vilà-Guerau de Arellano et al., 2004). Equation (8) shows a clear relationship between ER_atmos and ER_forest following the mixed-layer theory.

Using the definition of Eq. (1), we can extend Eq. (8) to include the effect of advection of O₂ (adv${}_{{\mathrm{O}}_{\mathrm{2}}}$) and CO₂ (adv${}_{{\mathrm{CO}}_{\mathrm{2}}}$), which is, next to entrainment, the second important large-scale process influencing the O₂ and CO₂ values:

$$\begin{array}{}\text{(9)}& {\mathrm{ER}}_{\mathrm{atmos}}={\mathrm{ER}}_{\mathrm{forest}}\cdot \left({\displaystyle \frac{\mathrm{1}+{\mathit{\beta }}_{{\mathrm{O}}_{\mathrm{2}}}+\frac{h}{(F_{{\mathrm{O}}_{\mathrm{2}}}{)}_{\mathrm{s}}}\cdot {\mathrm{adv}}_{{\mathrm{O}}_{\mathrm{2}}}}{\mathrm{1}+{\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}+\frac{h}{(F_{{\mathrm{CO}}_{\mathrm{2}}}{)}_{\mathrm{s}}}\cdot {\mathrm{adv}}_{{\mathrm{CO}}_{\mathrm{2}}}}}\right).\end{array}$$

Note that in this paper, we mostly focus on cases without advection. We include it here for completeness and discuss the influence of advection in Sect. 5.2.

In Appendix A1 we analyse Eq. (8) by determining when ER_atmos would theoretically be close to ER_forest during the day. We show that the β values are of particular importance here: when the β's of O₂ and CO₂ are equal or very small, ER_atmos gives the same signal as ER_forest. To fully unravel the diurnal variations of ER_atmos under realistic conditions and identify influencing factors, we need to analyse a real case. Therefore, we study two observed situations by means of the CLASS coupled land–atmosphere model, which we will describe in Sect. 3.2.

3 Methods

In this section we describe the measurements that were used in this study, together with the mixed-layer model used to evaluate the ER_atmos and ER_forest signals.

3.1 Hyytiälä 2018 and 2019 measurement campaigns

The observational data were obtained from the SMEAR II Forestry Station of the University of Helsinki in Finland, located in Hyytiälä, Finland (61°51^′ N, 24°17^′ E, +181 MSL) (Hari et al., 2013). The SMEAR II station serves as a measurement site within a boreal forest equipped with a 128 m tower for continuous measurements of atmospheric variables, fluxes, and greenhouse gas mole fractions. These data are accessible at https://smear.avaa.csc.fi/ (last access: 24 June 2024). The tower is situated in a homogeneous Scots pine forest, with an average canopy height of 18 m and podzolic soil. The measurement site is predominantly influenced by the surrounding forest and is minimally impacted by signals of fossil fuel combustion (Faassen et al., 2023). For a comprehensive description, see Hari et al. (2013).

During the spring and summer of 2018 (3 June until 2 August) and 2019 (10 June until 17 July), two measurement campaigns, referred to as OXHYYGEN (Oxygen in Hyytiälä), were conducted at Hyytiälä. Continuous measurements of both O₂ and CO₂ mole fractions were taken at two heights (125 and 23 m). O₂ was measured using an Oxzilla II fuel cell analyser, and CO₂ was measured with a non-dispersive infrared (NDIR) photometer (URAS26). Further details on these measurements and the measurement systems are given in Faassen et al. (2023). The measurement precision for O₂ was 19 per meg, and for CO₂ it was 0.07 ppm. Although the precision for O₂ is relatively poor compared to previous studies, it is still adequate for studying the diurnal timescale, as shown in Faassen et al. (2023).

O₂ measurements are typically expressed as δO₂ / N₂ ratios in per meg units due to the high abundance of O₂ in the atmosphere (20.946 %) classifying it as a non-trace gas. For direct comparison with CO₂ and implementation into our model, we convert per meg to ppm equivalents (ppmEq) by multiplying with the standard mole fraction of O₂ in air of 0.20946 (Keeling et al., 1998).

During the OXHYYGEN campaigns, radiosondes were launched on multiple days several times per day to quantify the impact of boundary layer dynamics on the O₂ and CO₂ diurnal cycles. The radiosondes (Windsond, model S1H3-R, Sweden) measured vertical profiles of air pressure, wind speed, wind direction, relative humidity, and temperature, with flight heights reaching a maximum of 4500 m and rising rate of about 1.7 m s⁻¹. The measurements have an accuracy of 1.0 hPa for air pressure, 5 % for wind speed, 0.2 °C for temperature, and 1.8 % for relative humidity. The temperature and humidity probe has a response time of 6 s, effectively averaging over about 10 m of altitude. For our analysis, we computed vertical profiles of potential temperature (θ) and specific humidity (q) based on pressure, temperature, and relative humidity measurements. Based on the vertical profile of vertical temperature, we also determine the boundary layer height with the parcel method (Kaimal and Finnigan, 1994). Figure 2 shows examples of vertical profile measurements of θ for 12 July 2019 and 24 July 2018.

3.2 Modelling setup in CLASS

3.2.1 Implementation of CO₂ in CLASS

CLASS serves as a fundamental tool that enables further understanding of specific processes within the atmospheric boundary layer. Several studies have shown that CLASS is successful in reproducing observational data (Vilà-Guerau de Arellano et al., 2012, 2019; Schulte et al., 2021). The study of Ouwersloot et al. (2012) specifically showed that CLASS is able to reproduce the boundary dynamics at the Hyytiälä measurement site. Within CLASS, the vegetation is described using a big-leaf model. The surface stomatal conductance that is representative for the canopy is up-scaled from leaf stomatal conductance by integrating over the leaf area index and incorporating soil moisture. The leaf stomatal conductance is calculated with the A-g_s model. The A-g_s model relates leaf stomatal conductance (g_s) to the net leaf CO₂ assimilation (A) (Jacobs et al., 1996; Ronda et al., 2001). The model computes the dependence of g_s and A on the internal CO₂ mole fraction, the amount of light, the atmospheric temperature, the vapour pressure deficit, and the soil water content at the root zone. Finally, the canopy net CO₂ assimilation is obtained with a function that is inspired by Fick's law of diffusion, based on the difference in the atmospheric CO₂ and internal CO₂ mole fractions, the aerodynamic resistance, and the surface stomatal conductance. The soil respiration is implemented as a function of soil temperature and soil moisture (Vilà-Guerau de Arellano et al., 2012). Combining the net assimilation (A_n) of the plants at canopy level and the soil respiration flux results in the net ecosystem exchange (NEE). This means that the model does not produce exactly the GPP and TER fluxes. The differences between A_n and GPP, and soil respiration and TER, are not directly relevant for our study and we therefore refer to GPP and TER in the following sections, as these terms are more commonly used in the atmospheric CO₂ community. The water cycle is connected to the CO₂ cycle through the surface stomata and the soil moisture inhibition functions for assimilation and respiration.

3.2.2 Implementation of O₂ in CLASS

To model both ER_forest and ER_atmos, we incorporated the surface flux and the atmospheric mole fraction of O₂ into the CLASS model. We represent the surface flux of O₂ by multiplying the ER of assimilation (ER_a) and the ER of respiration (ER_r) with the CLASS-calculated CO₂ fluxes at the canopy scale. We used the observationally derived ER_a and ER_r values as previously determined in Faassen et al. (2023) for the same site, which were 0.96 and 1.03, respectively. The net surface flux of O₂ was then resolved with the following equation:

$$\begin{array}{}\text{(10)}& F({\mathrm{O}}_{\mathrm{2}}{)}_{\mathrm{s}}=F_{{\mathrm{CO}}_{\mathrm{2}}\left(\mathrm{a}\right)}\cdot -{\mathrm{ER}}_{\mathrm{a}}+F_{{\mathrm{CO}}_{\mathrm{2}}\left(\mathrm{r}\right)}\cdot -{\mathrm{ER}}_{\mathrm{r}},\end{array}$$

where F(O₂)_s is the net O₂ surface flux above the canopy, $F_{{\mathrm{CO}}_{\mathrm{2}}\left(\mathrm{a}\right)}$ is the net assimilation flux, and $F_{{\mathrm{CO}}_{\mathrm{2}}\left(\mathrm{r}\right)}$ is the soil respiration flux. The change of atmospheric O₂ over time was resolved with Eq. (A1) (similar to Eq. 1) and the entrainment flux is based on Eq. (A2) (see also Eq. 2). Note that the ER_a from Faassen et al. (2023) was based on GPP fluxes, and this ER_a is now linked to the net assimilation flux (GPP minus the photo- and dark respiration) of the model (Jacobs et al., 1996; Ronda et al., 2001). Seibt et al. (2004) and Ishidoya et al. (2013) showed that ER_a values based on net assimilation have similar values compared to the 0.96 based on GPP. We therefore do not expect this discrepancy to influence our results.

Figure 3 Schematic overview of how two processes with different ER signals produce a combined ER signal that is not necessarily the average of the two processes nor necessarily falls within the range of the two combined ER signals. This is due to the different signs of the O₂ and CO₂ fluxes. The example is given for combining the ER signal of assimilation (ER_a) and respiration (ER_r) into ER_forest and uses values from our study that are by coincidence larger and smaller than 1.

It is important to note that the resulting ER_forest signal is not the (weighted) average between ER_a and ER_r, as was also shown by Faassen et al. (2023). The ER_forest signal results from the TER and GPP fluxes with different sizes and signs, each with their own ER signals (ER_r and ER_a, respectively). Figure 3 shows that the resulting ER_forest signal does not necessarily fall within the range of the ER_a and ER_r signals, because the TER and GPP have opposite signs of the O₂ and CO₂ fluxes. This counter-intuitive situation can also occur for combining signals with different isotopic signatures (Miller and Tans, 2003).

3.2.3 Initial conditions

We determined initial and boundary conditions for two cases to constrain the model to the observations. One case was based on the year 2019 (base case) and the other case was based on the year 2018 (characterized by a warm summer in Finland; Peters et al., 2020; Lindroth et al., 2020). Using the two years to initialize CLASS, we were able to better constrain the vegetation's response in the CLASS model under extreme conditions. For each year, we selected one representative day for initialization and validation of the CLASS model. We used 10 July 2019 for the base case and an aggregate between 28 and 29 August 2018 for the warm case. The initial and boundary conditions for initialization of the CLASS runs can be found in Tables C2 and C3 in the Appendix. Note that the initial jumps (Δ_(ft-bl)) of O₂ and CO₂ are based on the best fit between the model and the observations during the day, as direct observations of the jumps were not available. A detailed discussion can be found in Sect. 5.3.

We deliberately made only minimal adjustments for the initialization of the 2018 case compared to the 2019 base case to ensure consistency. We assumed that the initial relative humidity remained constant at 80 % regardless of temperature variations, similar to the studies of Vilà-Guerau de Arellano et al. (2012) and van Heerwaarden and Teuling (2014).

We adjusted several parameters of the A-g_s land–surface scheme and the soil respiration to improve the agreement between the surface fluxes of the model and the observations in Hyytiälä for both the base case (2019) and the warmer case (2018) (Table C2). We decreased the mesophyll conductance (g_m: 2 mm s⁻¹) to better match pine forest conditions (Gibelin et al., 2008; ECMWF IV, 2014; Visser et al., 2021). Furthermore, the reference temperature of g_m ($T_{\mathrm{2}\left(g_{\mathrm{m}}\right)}$: 305 K) was increased to reduce afternoon plant stress and to make the CLASS run more comparable with the observations. Lastly, we adjusted the curvature of the drought-response curve (c_β) from zero to 15 % (Combe et al., 2015), given that several studies have demonstrated the pine forest in Hyytiälä to be relatively resilient to lower soil moisture values, thus necessitating a higher (c_β) value (Gao et al., 2017; Lindroth et al., 2020).

3.2.4 Sensitivity analyses

We conducted two sensitivity analyses to gain a deeper understanding of the ER_atmos behaviour under varying conditions and to identify factors that lead to a smaller difference between ER_atmos and ER_forest. Specifically, we looked at changes in ER_atmos resulting from changing the different components of Eq. (8). The first sensitivity analysis uses the 2019 base case and investigates the effect of background air with a different composition by altering the initial jumps of O₂ and CO₂. By changing only the initial jump and keeping the rest of the 2019 case the same, we simulate situations in which the free troposphere mole fractions of O₂ and CO₂ have changed. In the second sensitivity analysis we examined the impact of climate conditions by modifying the soil moisture and air temperature, mimicking the conditions observed during the 2018 heatwave. Table C1 presents the variables used for initializing four cases for these two sensitivity studies.

4 Results

In this section, we first show our results for the validation of the CLASS model using observations (Sect. 4.1). Subsequently, we discuss the diurnal variability of both the ER_forest and ER_atmos signals (Sect. 4.2). We identify three distinct periods throughout the day in which ER_atmos shows large variability (Sect. 4.2.1). We address the large ER_atmos that we find in both the observations and the model results (Sect. 4.2.2). Finally, we perform sensitivity analyses to study the effects of changing large-scale conditions, in order to show that our findings are not only valid for a single day (Sect. 4.3).

4.1 Validation of the O₂ and CO₂ model results

Overall, the modelled O₂ and CO₂ diurnal cycles match well with the observational data. Figures A3 and A2 in Appendix A3 show that CLASS accurately reproduces the diurnal cycles and captures the O₂ mole fraction changes on a daily timescale for both 2018 and 2019 (Fig. A3b and c). The figure shows that the differences between the two years are relatively small and indicate that the boundary layer dynamics and the surface fluxes are well represented in CLASS. To accurately replicate the rapid decrease of CO₂ and the sharp increase of O₂ during the rapid growth of the atmospheric boundary layer (between 06:30 and 11:30), we adjusted the jump between the boundary layer and the free troposphere (Δ_(ft-bl)) for both O₂ (30 ppmEq) and CO₂ (8 ppm), ensuring that the model aligned with the measurements. Based on values from previous studies, it is realistic for the CO₂ jump to range between 8 and 40 ppm (Vilà-Guerau de Arellano et al., 2004; Casso-Torralba et al., 2008). While there are limited data available to validate the jump of O₂, based on preliminary results of a campaign in Loobos, the Netherlands, a jump of 30 ppmEq for O₂ seems reasonable. Our chosen combination of O₂ and CO₂ jumps remains an uncertain component in our analysis and will be further discussed in Sect. 5.3.

4.2 Diurnal variability of ER_atmos and ER_forest in 2018 and 2019

In this section we discuss the diurnal variability of the ER_atmos signal for the 2018 and 2019 cases. First we focus on the budget components (GPP, TER, and entrainment) that influence the O₂ and CO₂ signals (Sect. 4.2.1). To complete the analysis, we support the numerical analysis with Eq. (8) to gain a more comprehensive understanding of the underlying processes driving the ER_atmos signal for the 2019 case (Sect. 4.2.2).

Figure 4 Diurnal cycles of O₂ and CO₂ mole fractions (a) and ER_atmos and ER_forest (b) as modelled with CLASS for the selected days in 2018 and 2019. We identify three distinct periods based on panels (c) and (d), which show the tendencies for the 2019 case (change over time) for CO₂ and O₂ for each process that influences their mole fractions (Eq. 1): P1 05:00–06:30, P2 06:30–11:30, and P3 11:30–18:30 LT. The symbols represent half-hourly averaged values of the CLASS model output.

4.2.1 The three distinct periods of the ER_atmos signal during daytime

The ER_atmos signals obtained for the 2018 and 2019 experiments display large variability throughout the daytime (panels a and b in Fig. 4). We identify three distinct periods during the day based on the processes shown in Fig. 4c and d: (1) the early morning regime (P1, 05:00–06:30 LT) characterized by an increasing net CO₂ flux into the forest but a non-growing boundary layer (Fig. A3a), during which the ER_atmos signal during P1 is still relatively close to ER_forest; (2) the entrainment-dominant period (P2, 06:30–11:30 LT), where air from a residual layer or air masses from the free troposphere are entrained into the boundary layer and significantly influence the signals, leading to large ER_atmos values with an average greater than 3 and extreme values reaching close to 5; (3) the afternoon period (P3, 11:30–18:30 LT), where surface processes dominate the observed signals and ER_atmos moves slowly back towards ER_forest and becomes more consistent with values expected for surface processes. The ER_atmos values during the three identified periods show good agreement between the observations and the model results (Table 1). This analysis confirms from a model perspective that values above 2 for ER_atmos, as we reported in Faassen et al. (2023), are indeed possible. Figure 4c and d give first indications of what could cause these high values for ER_atmos: high influence of entrainment and a different behaviour of the tendencies that influence O₂ compared to CO₂. In the next section we discuss the diurnal behaviour of ER_atmos in more detail by using Eq. (8).

We find that ER_forest is much less variable throughout the day than ER_atmos (Fig. 4b). In the early morning and later afternoon, the ER_forest value is lower than during the midday period. This is caused by a TER flux (with a higher ER signal) almost equal to the GPP flux (with a lower ER signal) caused by low sunlight (Fig. 3). At midday, the assimilation of CO₂ by the canopy, with a lower ER signal, becomes increasingly dominant, causing the ER_forest signal to move closer to the ER_a value.

Table 1 ER_atmos (calculated as the slope of the O₂ and CO₂ mole fractions) and ER_forest values for the selected days in 2018 and 2019 for both observations (Obs) and the CLASS model for the three selected periods (P1: 05:00–06:30, P2: 06:30–11:30, and P3: 11:30–19:30 LT). The uncertainties in the observed ER_atmos and ER_forest signals are determined following Faassen et al. (2023). Note that due to limited observational data, we were unable to derive ER_atmos values for P1 and P2 in 2018 and for P1 in 2019. NA signifies not available.

Figure 5 Diurnal variability of the different components of Eq. (8) for the base case (2019) and the warm case (2018) derived using the CLASS model. Panels (a) and (e) show the β values for CO₂ and O₂, where β is the entrainment flux divided by the surface flux (Eq. 8); panels (b) and (f) show the net surface flux; panels (c) and (g) show the jumps between the free troposphere and the boundary layer (Δ_(ft-bl)); and panel (d) shows the entrainment velocity (w_e). The vertical lines represent three distinct periods: 05:00–06:30 (P1), 06:30–11:30 (P2), and 11:30–18:30 LT (P3).

4.2.2 Explanation of the large ER_atmos values

Analysing the diurnal cycle of the different components of Eq. (8) for the 2019 case reveals that the peak value of ER_atmos during P2 is caused by the higher β values (the entrainment flux divided by the surface flux) for O₂ compared to CO₂ (Fig. 5). The difference between ${\mathit{\beta }}_{{\mathrm{O}}_{\mathrm{2}}}$ and ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$ is a result of the high Δ_(ft-bl)O₂ $/$ Δ_(ft-bl)CO₂ ratio (higher than 3). The terms Δ_(ft-bl)O₂ and Δ_(ft-bl)CO₂ represent the jump across the boundary layer top, and each has a different diurnal cycle caused by a different surface flux (Fig. 5c and g). The different diurnal cycles for the jumps lead to an increase in the Δ_(ft-bl)O₂ $/$ Δ_(ft-bl)CO₂ ratio, consequently raising the ratio between the β values. This effect is further amplified by a higher surface flux of CO₂ compared to O₂, caused by an ER_forest value that is slightly lower than 1. The peak value of ER_atmos during P2 occurs when both w_e and the Δ_(ft-bl)O₂ $/$ Δ_(ft-bl)CO₂ ratio are high and the surface fluxes are still relatively low. This combination contributes to the distinctive peak in ER_atmos observed during P2.

Later in the afternoon (P3), both β values gradually decrease and become similar, resulting in an ER_atmos signal that becomes closer to ER_forest. This indicates that ER_atmos becomes more representative for surface processes (see also Appendix Sect. A1). This decrease in P3 is primarily caused by a reduction in the entrainment velocity (w_e) (Fig. 5d), indicating slow growth of the atmospheric boundary layer at end of the day (Fig. A3). Additionally, the β values become more similar because Δ_(ft-bl)O₂ moves closer to Δ_(ft-bl)CO₂ during this period (Fig. 5c and g), as caused by the mixing of air with the surface.

The ER_atmos signals exhibit higher values than the theoretical analysis in Appendix Sect. A1, because the diurnal cycles of the components of Eq. 8 are taken into account (Fig. A1 vs. Fig. 5). Each component of Eq. (8) follows its individual diurnal cycle, leading to higher ER_atmos values. Consequently, ER_atmos integrates individual contributions of several processes, particularly during P2, since it is dominated by the influence of mixing with large-scale processes. Careful consideration is needed when interpreting the ER_atmos signal during this period. During P3 the ER_atmos signal appears to align with ER_forest at the end of the day. However, in the 2019 case, this alignment was only observed for a very short period.

We find only small differences in the diurnal behaviour of the ER_forest and ER_atmos signals between the 2018 and 2019 cases (Figs. 4 and 5). The ER_forest value is lower in 2018 compared to 2019, specifically in the early morning and at the end of the day. This can be attributed to a higher respiration flux caused by the elevated air and soil temperatures during that day in 2018 (Fig. A3e). A higher TER flux compared to the GPP flux will decrease the ER_forest value (Fig. 3). While we do not have direct measurements of ER_r and ER_a for 2018 and 2019, it is likely that the overall diurnal cycle pattern of ER_forest in Fig. 4b (low ER_forest values in the morning and afternoon, higher ER_forest values during midday) for both years would have remained consistent. Previous studies suggest that ER_r is generally higher than ER_a, even under different atmospheric conditions (Angert et al., 2015; Fischer et al., 2015; Hilman et al., 2022). The effect of a warmer and dryer environment on the ER_atmos signal will be further quantified in Sect. 4.3.2 with a more extreme case.

4.3 Sensitivity analyses: effects of changing large-scale conditions

With the next two sensitivity analyses, we evaluate whether our findings for the 2019 case are exceptional or whether they can also occur under different (large-scale) conditions. We therefore analyse days with different initial conditions compared to our 2019 and 2018 cases. We focus on the effect of changes in background air (Sect. 4.3.1) and the effect of changes in climate conditions (soil moisture and air temperature, Sect. 4.3.2). With these sensitivity analyses, we show the complexity of the ER_atmos signal and all the processes by which it can be influenced. Fig. 6 is used to illustrate how ER_atmos is formed by the different components of Eq. (8).

Figure 6 Components of Eq. (8) and how these influence the ER_atmos signal, including the exchange ratio of the forest (ER_forest), the ratio between the net surface flux (F_s) and the entrainment flux (F_entr) which result in β, the jump between the free troposphere and the boundary layer (Δ(ft-bl)), and the entrainment velocity (w_e). The right part of the figure shows the variables that are changed in the two sensitivity analyses: the background air in the free troposphere ([O₂](ft) and [CO₂](ft)) and the initial soil moisture index (SMI) in combination with a high initial potential temperature (θ₀), which will influence the ratio between the sensible heat flux (SH) and the latent heat flux (LH) at the surface. The dotted arrows indicate a negative influence and the solid arrows indicate a positive influence.

4.3.1 Effects of changing background air on ER_atmos

Changing the background air in the free troposphere by decreasing the initial jump ratio or the jump sizes of O₂ and CO₂ compared to the 2019 case moves the ER_atmos signal closer to ER_forest during P2 and P3 (Fig. B1). A lower jump ratio than in the 2019 case but still relatively high jump values (Δ_(ft-bl)O₂=30 ppmEq and Δ_(ft-bl)CO${}_{\mathrm{2}}=-\mathrm{20}$ ppm) lead to a decrease in the peak of ER_atmos during P2 and bring ER_atmos closer to ER_forest during P3 (yellow line in Fig. B1). As the jump ratio decreases, ${\mathit{\beta }}_{{\mathrm{O}}_{\mathrm{2}}}$ becomes less dominant and closer to ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$. When the O₂ and CO₂ β values become closer, the ER_atmos value also moves closer to ER_forest (Fig. 6). However, this does not necessarily mean that the surface has become more dominant, since the Δ_(ft-bl) values are still relatively high.

Reducing the jump sizes of both O₂ and CO₂ (Δ_(ft-bl)O₂=10 and Δ_(ft-bl)CO${}_{\mathrm{2}}=-\mathrm{8}$) still results in a relatively high peak for ER_atmos during P2 and brings ER_atmos closer to ER_forest during P3 (purple line in Fig. B1). Including the diurnal cycle of the jumps accounts for the effect that the CO₂ jump changes from a negative to a positive value during the day. When the initial CO₂ jump is lower, the sign change occurs earlier in the day and leads to a more negative ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$ value. This leads to higher ER_atmos values during P2 (Fig. 6). In contrast, a lower jump size would cause the ER_atmos signal to move more quickly towards ER_forest during P3, because the surface fluxes dominate over the lowered entrainment flux.

Guided by our theoretical and numerical results and constrained by observations, a high ER_atmos signal during the entrainment-dominant period (P2) can therefore be a result of two cases:

1. The Δ_(ft-bl)O₂ is substantially larger compared to Δ_(ft-bl)CO₂ and, therefore, ${\mathit{\beta }}_{{\mathrm{O}}_{\mathrm{2}}}$ dominates over ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$.

2. Δ_(ft-bl)CO₂ changes sign from negative to positive and, as a result, ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$ becomes negative, resulting in a denominator closer to zero.

Changes in the background air result in a distinct change in the diurnal pattern of ER_atmos. The difference between the ER_atmos and ER_forest signals could therefore provide extra information on the changes of large-scale processes. This is further discussed in Sect. 5.2.

4.3.2 Effect of climate conditions on ER_atmos and ER_forest

By studying the influence of changes in air temperature and soil moisture index (SMI: [soil moisture – w_wilt]/[w_fc−w_wilt]) on the ER_atmos signal (see Fig. 7), we gain insights into how climate conditions can effect ER_atmos compared to ER_forest. This allows us to study the effects of seasonality or future climate with dryer and warmer conditions. The 2018 case already showed how the ER_atmos signal could change with decreasing SMI and increasing temperature compared to a more normal year in 2019 (Figs. 4 and 5). As a next step, we evaluate the full range of how ER_atmos could change and how ER_atmos compares to ER_forest. Given the same net radiation, a higher SMI enhances soil respiration, photosynthesis, and latent heat fluxes, and thus decreases the sensible heat flux because of the energy balance closure. This therefore leads to smaller boundary layer growth and, as a result, decreases the entrainment velocity (see Fig. 6). In addition, higher air temperatures accelerate both photosynthesis and respiration up to a threshold (Jacobs et al., 1996), resulting in increased GPP and TER fluxes. A lower SMI in combination with higher temperatures can stress plants, leading to decreased O₂ and CO₂ surface fluxes and an enhanced sensible heat flux. This will increase boundary layer growth and the entrainment velocity (Eqs. 2 and 4). Note that there are also minor changes in ER_forest when the SMI and air temperature change as a result of GPP and TER changes.

Figure 7 Evaluation of the ratio between ER_atmos and ER_forest as a function of two key variables that show the effect of a drier and warmer climate: the soil moisture index (SMI) and the initial potential temperature (θ₀). Two periods in the day are analysed: (a) during the maximum value of w_e at 08:14 LT (P2) and (b) at the end of the day between 15:30 and 17:30 LT, when w_e is minimal (P3). The black lines in panel (a) indicate ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$, which is the ratio between the entrainment and the surface flux. The black lines in panel (b) indicate net CO₂ surface flux values in ppm m s⁻¹. The coloured symbols (brown and light blue) indicate the example cases that are also shown in Fig. B2, the black dot is the 2018 case, and the dark blue dot the 2019 case.

Increasing or decreasing the SMI in combination with changes in air temperature makes the diurnal variability of ER_atmos more complex, because all components of Eq. (8) are now affected (Figs. 6 and 7). We focus on two particular locations in the parameter space shown in Fig. 7: a low soil moisture (red symbol) and a high soil moisture case (green symbol), both with higher temperatures compared to the 2019 case (Fig. B2).

A lower soil moisture of 0.14 m³ m⁻³ (SMI =0.27) with an air temperature of 290 K decreases ER_atmos during P2 and increases ER_atmos during P3 compared to the 2019 base case (red lines in Fig. B2 and red symbol in Fig. 7). The lower ER_atmos values during P2 are primarily a consequence of a more dominant entrainment flux. Due to a decrease in the O₂ and CO₂ surface fluxes because of stressed plants, both the Δ_(ft-bl) values for O₂ and CO₂ change more slowly and remain high. Higher Δ_(ft-bl) values, along with a higher entrainment velocity caused by higher sensible heat flux, lead to elevated entrainment fluxes. By increasing both the O₂ and CO₂ entrainment fluxes and decreasing the O₂ and CO₂ net surface fluxes, the β values increase and the ratios of the β values move towards the Δ_(ft-bl) ratios. As a result, ER_atmos also moves towards the Δ_(ft-bl) ratios multiplied with the ER_forest signal (Fig. 6). This is similar to the effect observed when increasing the initial jumps of both O₂ and CO₂ (Sect. 4.3.1). The β values stay high during P3 because of the low net O₂ and CO₂ surface fluxes. Therefore, the ER_atmos signal also remains close to the ratio of the Δ_(ft-bl) values during P3 and the ER_atmos signal does not approach ER_forest (Fig. 6).

In contrast, a higher soil moisture of 0.22 m³ m⁻³ (SMI =0.64) with an air temperature of 290 K increases the ER_atmos signal during P2 and decreases the ER_atmos signal during P3 compared to the 2019 base case (green lines in Fig. B2 and green symbol in Fig. 7). This is consistent with the effect observed when lowering the initial Δ_(ft-bl) value (Sect. 4.3.1).

In addition to the conclusions in Sect. 4.3.1 on the causes of the high ER_atmos signals during P2, the sensitivity analyses for changing climate conditions showed that the large differences between ER_atmos and ER_forest at the end of the day (P3) can be caused by

1. A substantially larger Δ_(ft-bl)O₂ compared to Δ_(ft-bl)CO₂, causing ${\mathit{\beta }}_{{\mathrm{O}}_{\mathrm{2}}}$ to dominate over ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$.

2. High ${\mathit{\beta }}_{{\mathrm{O}}_{\mathrm{2}}}$ and ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$ values because of high O₂ and CO₂ entrainment fluxes and/or low net O₂ and CO₂ surface fluxes.

Our two sensitivity analyses show that several factors, including the entrainment velocity, the Δ_(ft-bl) values and their ratio, and the net surface flux of CO₂, can significantly influence the diurnal behaviour of ER_atmos. When using ER_atmos as an indication of ER_forest, these four factors should be carefully considered. This is crucial to correctly interpret ER_atmos values and to understand the underlying processes that influence the carbon exchange above a forest canopy.

5 Discussion

In this section we first address the evaluation of the CLASS model (Sect. 5.1). Secondly, we elaborate on the issues we found with ER_atmos and how this value should (and should not) be used (Sect. 5.2). Thirdly, we discuss the importance of the differences between the free troposphere and boundary layer values for O₂ and CO₂ (Sect. 5.3). Finally, we put our work in perspective by comparing it to other studies using atmospheric O₂ (Sect. 5.4) and to studies on other carbon cycle tracers (Sect. 5.5).

5.1 Evaluation of the CLASS model

Our implementation of O₂ in the CLASS model could be improved in future studies. Similar to the approach used by Yan et al. (2023), both the ER_r and ER_a signals were kept constant and did not account for potential variations under different climate conditions. To advance our understanding of the ER signals over forest canopies, it is crucial to incorporate ER signals that can respond to varying soil and atmospheric conditions. For instance, the ER_r of soil respiration depends on air temperature and soil moisture (Hilman et al., 2022; Angert et al., 2015), while ER_a is primarily influenced by light at leaf level and nitrogen availability in the soil (Bloom, 2015; Fischer et al., 2015). Additionally, in our current implementation, we did not include the ER for stem respiration (ER_stem) (Hilman and Angert, 2016) due to the absence of stem respiration in the CLASS model.

While we utilized CLASS in this study as a proof of concept to demonstrate how ER_atmos can change during the day; employing a more elaborate model could allow for more detailed exploration of these ER_atmos dynamics and the contributions of various processes. Models with more vertical levels could simulate vertical gradients and analyse differences in the ER_atmos signal at various heights, similar to the approach used in Yan et al. (2023). Implementing more vertical levels provides the opportunity to determine the dominance of large-scale processes over small-scale surface processes at different measurement heights. By incorporating a canopy into the model, the surface resistance could be accounted for, enhancing the accuracy of the modelled surface fluxes. Furthermore, exploring larger temporal and spatial scales could yield valuable insights into the variability of ER_forest over time and space, in contrast to our CLASS model that is only valid during the day when the SH flux is larger than zero. Increasing the temporal scale gives the opportunity to improve estimates of ER_forest. This also has the potential to improve estimates of the global biospheric ER, currently taken to be 1.1 (Severinghaus, 1995).

5.2 How ER_atmos should be used

Single-height O₂ and CO₂ measurements and their ER_atmos signal should be analysed very carefully when using ER_atmos as an indicator of surface exchange. During the complete diurnal cycle, ER_forest should be utilized as the primary indicator of the ER signals from the surface, while ER_atmos should not be used for this purpose. In situations where only one height measurement is available, and therefore only ER_atmos can be obtained, a first estimate of ER_forest could be made using ER_atmos. The ER_atmos signal at the end of the day should then be used to avoid the large influence of entrainment earlier in the day. However, any analysis or discussion based on this estimation should include a comprehensive examination of how entrainment might have influenced the ER_atmos signal. This also applies for less representative or non-typical days where the mixed-layer theory may be difficult to apply. An example of such a case is given by Casso-Torralba et al. (2008), where it is shown that entrainment is still important on a non-typical day when polluted air influences the diurnal CO₂ measurements.

Several studies have shown that ER_atmos can also serve as an indicator of potential advection from carbon source/sink regions (Ishidoya et al., 2020, 2024). However, caution should be exercised when directly inferring the specific source based solely on the ER_atmos value. Equation (9) shows that mixing advected air with the air above a forest will result in an ER_atmos signal that cannot be directly linked to the source of the advected air. This is because mixing two ER signals with opposite fluxes does not result in a weighted average (Fig. 3). Advection of a source with a known ER signal but with different magnitudes can therefore give different ER_atmos values. A solution could be to include other tracers in the analysis such as NO_x or CO (Liu et al., 2023a).

When two or more measurement heights of O₂ and CO₂ are available and ER_forest can therefore be derived, the ER_atmos of a single height could be used to provide extra information on large-scale processes by analysing the difference between ER_atmos and the ER_forest signal. During the day, ER_atmos provides insights into larger-scale processes, while ER_forest reflects local or small-scale processes. Therefore, any discrepancy between ER_atmos and ER_forest indicates a significant influence of large-scale processes. Nonetheless, the exact difference between ER_atmos and ER_forest should not be used as an indication of the strength of the influence of large processes. To get more detail on how the large scale processes change between days, the diurnal cycle of ER_atmos has to be compared during the entrainment dominant period (P2) and the surface-dominant periods (P3). During P2, an increase in the difference between ER_atmos and ER_forest may be due to either a low ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$ or a change in the jump (Δ_(ft-bl)) ratio. If the cause is the former (low ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$), the ER_atmos signal during P3 should be closer to ER_forest. If the latter applies (a high jump ratio), ER_atmos should remain well above ER_forest in P3.

Figure 8 Schematic overview of how different ratios of the jumps of O₂ (Δ_(ft-bl)O₂) and CO₂ (Δ_(ft-bl)CO₂) are formed in the nighttime and how the ratio relates to the exchange ratio of the forest (ER_forest). Panel (a) gives an explanation of the schematics and the other panels show four possibilities of different jump ratios (Δ_(ft-bl)O₂ $/$ Δ_(ft-bl)CO₂): (b) the jump ratio is equal to ER_forest, (c) the jump ratio is larger than ER_forest, (d) the jump ratio is much larger than ER_forest, (e) the jump ratio is smaller than ER_forest. The bold lines represent the vertical profile just after sunset and the shaded lines represent the vertical profile just before sunrise.

5.3 Different Δ_(ft-bl) ratios

Knowing the vertical profile of O₂ and CO₂, especially during sunrise, is essential to gain a more comprehensive understanding of the formation of different jump ratios (Δ_(ft-bl)O₂ $/$ Δ_(ft-bl)CO₂) and to better interpret the diurnal behaviour of the ER_atmos signal. However, due to a lack of observational data, we cannot validate the vertical profile of O₂ and CO₂ and the jump ratios. We therefore strongly recommend that future measurement campaigns include vertical measurements of both species. This can, for example, be achieved by flask sampling from aircraft, as we did in a recent campaign in the Netherlands, the preliminary results of which confirm that the values we have used here are realistic. Previous studies have also measured vertical profiles of O₂ and CO₂, but they primarily focused on well-mixed profiles during daytime or profiles over the ocean (Morgan et al., 2019; Stephens et al., 2021; Ishidoya et al., 2022). Hence, careful consideration of the timing and location of the vertical measurements is important to advance our knowledge of the diurnal behaviour of ER_atmos.

In the absence of observational data, we show with hypothetical situations that various jump ratios become possible (Fig. 8). Both the O₂ and CO₂ jumps are formed as a result of three processes: the mixed-layer value before sunset, the surface flux during the night, and the free troposphere value with the lapse rate (we assume the lapse rate to be 0 mol m⁻¹ for CO₂ and O₂). Most cases indicate that Δ_(ft-bl)O₂ is larger than Δ_(ft-bl)CO₂ above a forest, primarily because ER_forest is higher than 1.0 during the night (ER_r>1.0), as was shown in previous studies (Ishidoya et al., 2013; Angert et al., 2015; Hilman et al., 2022). It is noteworthy that the movement of the mixed-layer values from sunset to sunrise in Fig. 8 differs from its depiction in Fig. 2c, where the focus was primarily on the transition between sunrise and sunset. We ignore the effect of subsidence (caused by mesoscale or synoptic processes) on the jump evaluation in this analysis, because it is likely of less importance compared to the other three processes.

It is highly likely that the jump ratio between O₂ and CO₂ cannot be directly linked to a specific ER for a certain process because of the interplay between the three processes that form the O₂ and a CO₂ jump (Fig. 8d). The likelihood of both Δ_(ft-bl)O₂ and Δ_(ft-bl)CO₂ being zero at the end of the day is low, because the surface flux during the day would form a jump (Fig. 8c). Additionally, it is possible that Δ_(ft-bl)O₂ is smaller than Δ_(ft-bl)CO₂ at the end of the day (Fig. 8d) due to the daytime ER_forest being smaller than 1.0. Consequently, O₂ will exhibit a faster movement across the zero line, resulting in a significantly larger Δ_(ft-bl)O₂ $/$ Δ_(ft-bl)CO₂ ratio compared to ER_r.

Decoupling between the free troposphere and the boundary layer can lead to a scenario in which Δ_(ft-bl)CO₂ becomes larger than Δ_(ft-bl)O₂ (Fig. 8e). This can occur, for example, when the influence of fossil fuel sources causes a decrease in the O₂ mole fraction and an increase in the CO₂ mole fraction in the free troposphere, but large surface fluxes from the forest prevent such changes from occurring in the boundary layer. The jump ratio in this case again cannot be attributed to a single process. Some studies have demonstrated that decoupling between the boundary layer and the free troposphere can occur, leading to different ER signals (Sturm et al., 2005; van der Laan et al., 2014).

5.4 Comparison with other studies

To the best of our knowledge, no previous studies have reported such high deviations of ER_atmos from ER_forest or ER_atmos values higher than 2 for above-forest-canopy measurements, as we found in Faassen et al. (2023). Only Liu et al. (2023b) found a non-linear relationship between O₂ and other tracers that was difficult to explain. While some differences between ER_atmos and ER_forest have been observed in previous studies, these differences typically fall within a region of 0.5 (Seibt et al., 2004; Ishidoya et al., 2015; Battle et al., 2019; Yan et al., 2023). A possible reason for these smaller differences could be that most studies do not focus on such detailed diurnal analyses of ER_atmos for specific days, but rather aggregate data from multiple days, which could mitigate the extreme effects of entrainment by combining various jump possibilities. However, even in the study by Stephens et al. (2007), in which measurements at different heights are shown, no discernible difference in the ER_atmos signal for various diurnal cycles was observed, a finding that contrasts with our own analysis. The height at which measurements are made also influences the resulting ER_atmos signal. Closer to the canopy, the influence of entrainment is lower and ER_atmos is closer to ER_forest compared to measurements further away from the canopy (Faassen et al., 2023). However, we still found a high ER_atmos value of 2.28, even at a level just above the canopy (Faassen et al., 2023). Large values for ER_atmos have only been found at high-latitude measurement stations (Sturm et al., 2005), due to the influence of the ocean.

There are several possibilities that might explain a constant ER_atmos signal during the day which are not shown in our study. One possibility is that entrainment dominates throughout the day, caused by high jumps. If both the O₂ and CO₂ jumps are extremely high while the surface flux remains low, the ER_atmos value reflects the ratios between the jumps. In this scenario, ER_atmos cannot be used as an accurate indicator for the surface processes. Another explanation could be that the ER_forest signal is exactly 1.0 and entrainment is relatively low. When ER_forest equals 1.0, the diurnal cycle of the jumps would respond similarly. Together with a low entrainment flux (resulting from low jumps), it could lead to a constant ER_atmos value. Additionally, when the peak of ER_atmos occurs rapidly, there is a possibility that a low measurement precision would miss the extreme changes of ER_atmos. However, even in such cases, ER_atmos would still be influenced by entrainment, although its impact may be less discernible. It is crucial to note that in all these cases, ER_atmos remains influenced by entrainment to varying degrees.

Our study provides evidence that ER_atmos is almost always influenced by large-scale processes and their diurnal variability, specifically entrainment, making it important to exercise caution when using it as an indicator of surface ER processes. Instances where ER_atmos remains constant throughout the daytime and serves as a reliable indication for ER_forest are likely rare. In comparison to previous studies (Seibt et al., 2004; Stephens et al., 2007; Ishidoya et al., 2013; Battle et al., 2019), it is unclear why Faassen et al. (2023) yields such extreme values for ER_atmos while the other studies do not show this, even though our modelling study here confirms the extreme ER_atmos values. Therefore, we recommend conducting more studies or performing detailed analyses of existing O₂ and CO₂ datasets to gain a better understanding of how changes in ER_atmos vary with time and space.

5.5 Comparison with other multi-tracer analyses

The impact of changes in large-scale conditions such as entrainment on multi-tracer analyses above forest canopies extends beyond atmospheric O₂, encompassing other carbon cycle tracers such as carbon and oxygen isotopes (δ¹³C and δ¹⁸O) (Wehr et al., 2016) and carbonyl sulfide (COS) (Whelan et al., 2018). Caution is required when employing methods of determining ratios between two species (e.g. leaf relative uptake for COS and the ratios between different isotopes) that rely solely on single-height measurements. However, the influence of entrainment on these ratios would be less extreme compared to that on the ER_atmos signal, because both COS and isotopes move in the same direction as CO₂ itself. This is differs from the situation with O₂, which always moves in the opposite direction to CO₂. When both species that form the ratio move in the same direction, ratios of different processes could be averaged and a one-height measurement is more readily interpretable. Nevertheless, entrainment would still cause the two compounds that form the ratio to behave differently. We therefore emphasize the need to separately analyse the composition of the signal for each compound when ratios are analysed.

Furthermore, we demonstrate in this study the potential of using ER_atmos as an indicator of the extent of large-scale processes. Additional tracers can strengthen this approach. Also δ¹³C, δ¹⁸O, and COS signals exhibit differences between the surface and the free troposphere. Similar to O₂, the onset of entrainment causes these signals to mix, yielding insights into how large-scale processes influence the carbon cycle above a canopy (Berkelhammer et al., 2014; Vilà-Guerau de Arellano et al., 2019). By combining various tracers for CO₂, we can create a comprehensive picture of the effects of small- and large-scale processes that influence carbon exchange.

6 Conclusions

We used a mixed-layer model to analyse the diurnal behaviour of two exchange ratio (ER $=-{\mathrm{O}}_{\mathrm{2}}/{\mathrm{CO}}_{\mathrm{2}}$) signals above a forest canopy: the ER of the atmosphere (ER_atmos, determined from the change over time of O₂ and CO₂ mole fraction measurements at a single height above the canopy) and the ER of the forest (ER_forest, determined from O₂ and CO₂ fluxes derived from the vertical gradient observations at two levels). We disentangled the biophysical processes influencing ER_atmos to interpret single-height O₂ and CO₂ measurements and to evaluate how ER_atmos and ER_forest can be used to constrain carbon exchange above the canopy. The analysis is supported by the derivation of a new theoretical relationship that connects ER_atmos and ER_forest and by the use of a mixed-layer model that reproduces the O₂ and CO₂ diurnal cycles coupled to the dynamics of the atmospheric boundary layer. By combining the model with observations from a boreal forest during the two contrasting summers of 2018 and 2019, we found three distinct regimes during the day for ER_atmos.

We find that the entrainment of air from the free troposphere leads to a diurnal cycle in ER_atmos, resulting in three distinctive regimes: P1 at the start of the day, when the boundary layer has not yet started to grow; P2 when entrainment of air from the free troposphere into the boundary layer is dominant; and P3 at the end of the afternoon, when entrainment becomes negligible. ER_atmos can exhibit high values during P2 that cannot be attributed to an ER signal from a single process. During P3, ER_atmos becomes closer to ER_forest, and is therefore more representative of the forest exchange.

The large diurnal variability in ER_atmos shows that single-height O₂ and CO₂ measurements are insufficient as an indication of the O₂ $/$ CO₂ ratios of forest exchange. Our theoretical relationship between ER_atmos and ER_forest and the model results show that the large diurnal variability is a result of the different behaviours of the O₂ and CO₂ diurnal cycles, which results in ER_atmos values that cannot be attributed to a single process. To estimate the ER signal of the surface fluxes from above-canopy measurements, ER_forest should be used and, therefore, O₂ and CO₂ signals need to be measured at at least two heights to allow fluxes to be calculated from the vertical gradient. A single measurement height for O₂ and CO₂ could still be used to indicate the presence of advection of other carbon sources. However, the resulting ER_atmos signal should be analysed with care, taking into account the diurnal variability and the fact that the resulting ER is not necessarily the average of the individual ER signals of the contributing processes.

When O₂ and CO₂ measurements are available from two different heights, the relationship between ER_atmos and ER_forest during P2 and P3 could provide valuable information about the changes in large-scale carbon processes (e.g. entrainment) and their influence on the smaller-scale processes of the surface. A discrepancy between ER_atmos and ER_forest shows that large-scale processes occur together with small-scale processes at the surface. The difference between ER_atmos and ER_forest should be analysed with care, as the size of the difference is not a direct indication of the size of the influence of the large-scale processes. Differences between ER_forest and ER_atmos could be caused by several factors: changes in the size of the entrainment flux, the net surface flux or the difference between the free troposphere and the boundary layer (the “jump”) for O₂ and/or CO₂, or changes in the jump ratio between O₂ and CO₂.

In conclusion, single-height O₂ and CO₂ measurements need to be analysed with care, accounting for their dependence on canopy processes (represented by ER_forest) but also for their capacity to integrate large-scale processes, resulting in values that cannot be attributed to a single process. To represent the forest exchange, the ER_forest signal based on measurements at at least two heights should be used instead.

Appendix A: Appendix

A1 Evaluation of the theoretical relationship between ER_atmos and ER_forest

In this section, we analyse Eq. (8) to explore the response of ER_atmos to changes in the variables in this equation and to investigate when ER_atmos aligns with ER_forest and thereby accurately reflects local processes. Based on Eq. (8), the ER_atmos signal equals ER_forest when the β values of O₂ and CO₂ are equal. We can define four different regimes where the β values change significantly. As depicted in Fig. 1, we can define two regimes based on the entrainment velocity: an entrainment-driven (left panels in Fig. A1) and a photosynthesis-driven regime (right panels in Fig. A1). To complete the analysis, we considered two distinct cases for the jump of O₂ (top versus bottom panels).

Based on Eq. (8), we systematically varied Δ_(ft-bl)CO₂ and $(F_{{\mathrm{CO}}_{\mathrm{2}}}{)}_{\mathrm{s}}$ over plausible ranges and kept the other variables constant. As a result, we derived ER_forest $/$ ER_atmos ratios for these four regimes, where a value of 1.0 now indicates that ER_atmos is equal to ER_forest. The selected values and ranges for the four different cases were informed by initial conditions from the Hyytiälä case studied in Faassen et al. (2023) and the corresponding model simulations presented in Sect. 3.2.4.

Figure A1 Analysis of Eq. (8) for the entrainment- and photosynthesis-driven regimes. The ratio between ER_forest and ER_atmos is evaluated based on changes in ER_forest and the ratio of the jumps of O₂ and CO₂ between the free troposphere and the boundary layer (Δ_(ft-bl)) for four cases: with a high entrainment velocity (w_e=0.10 m s⁻¹) (a, c) and a low entrainment velocity (w_e=0.01 m s⁻¹) (b, d), and for situations with a high O₂ jump (Δ_(ft-bl)O₂=0.30 ppmEq) (a, b) and a low O₂ jump (Δ_(ft-bl)O₂=0.10 ppmEq) (c, d). The O₂ surface flux F(O₂)_s is kept constant for all the panels at 8.5 µmol m s⁻¹.

There are a few situations where the β values of O₂ and CO₂ are equal and these are indicated in Fig. A1 as the area between the black solid (ER_atmos deviates <1 % from ER_forest) and dashed lines (ER_atmos deviates <10 % from ER_forest):

1. During the photosynthesis dominant regime. When the entrainment velocity (w_e) is close to zero, both β values become zero. This is likely at the end of the day (right panels in Fig. 1).

2. When the β values for O₂ and CO₂ become equal, which happens when Δ_(ft-bl)O₂ $/$ Δ_(ft-bl)CO₂ = ER_forest. A specific case is when Δ_(ft-bl)O₂ = Δ_(ft-bl)CO₂. In this instance, ER_forest has to be 1.0 for the β values of O₂ and CO₂ to become equal. The β values of O₂ and CO₂ become closer during the lower O₂ jump case (lower panels in Fig. 1).

The last situation only occurs under very specific conditions when the ratio of the O₂ and CO₂ entrainment and surface fluxes are the same. This is visible in the left panels of Fig. A1, where only a small part of the graph shows values of ER_atmos close to ER_forest (indicated by the area between the solid lines). In contrast, during low entrainment velocities at the end of the afternoon, it is more likely that the ER_atmos values become close to ER_forest, and this is shown by the larger area in the right panels of Fig. A1. Low entrainment velocities could also occur when the growth of the boundary layer is reduced due to subsidence, although we do not focus on this specific case in this study.

There are also differences between ER_atmos and ER_forest that arise from variations in the β values. Figure A1 demonstrates that substantial differences between ER_atmos and ER_forest originate due to differences in the entrainment fluxes of both species. When Δ_(ft-bl)O₂ exceeds Δ_(ft-bl)CO₂, this implies a dominant entrainment flux of O₂ over CO₂ and ${\mathit{\beta }}_{{\mathrm{O}}_{\mathrm{2}}}$ deviates further from ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$ (Eq. 8). This effect is almost absent when the jumps themselves are lower, because the ER_atmos/ER_forest ratio stays around 1 (Fig.  A1c). Moreover, when Δ_(ft-bl)CO₂ transitions from negative to positive, the sign of ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$ also changes, subsequently elevating the ER_atmos values (Eq. 8).

Figure A2 Comparison between the 2019 and 2018 cases modelled with CLASS with the observational data for the potential temperature (θ) (a), specific humidity (q) (b), sensible heat flux (SH) (c), and latent heat flux (LH) (d).

Figure A3 Comparison between the 2019 and 2018 cases modelled with CLASS using the observational data for the boundary layer height (a), CO₂ (b), O₂ (c), the 2019 CO₂ surface fluxes (d), and the 2018 CO₂ surface fluxes (e).

ER_atmos can also become smaller than ER_forest when Δ_(ft-bl)CO₂ is larger than Δ_(ft-bl)O₂ (Fig. A1). This difference results in a large value for ${\mathit{\beta }}_{{\mathrm{CO}}_{\mathrm{2}}}$ compared to ${\mathit{\beta }}_{{\mathrm{O}}_{\mathrm{2}}}$, causing the ER_forest value to be multiplied by a factor of less than 1 and leading to a lower ER_atmos value than ER_forest (Eq. 8). By assessing ER_atmos and ER_forest values, we can see whether Δ_(ft-bl)O₂ exceeds Δ_(ft-bl)CO₂ (ER_atmos> ER_forest) or vice versa (ER_atmos< ER_forest).

This illustrative analysis, based on prescribed values in Eq. (8) and Fig. A1, provides an initial estimate of the variability in ER_atmos. However, it lacks insights into the diurnal behaviour of the individual components of Eq. (8) and their potential combinations.

A2 Implementation of O₂ in CLASS

The following equation shows the implementation of the tendency (change over time) of O₂ into CLASS:

$$\begin{array}{}\text{(A1)}& {\displaystyle \frac{\mathrm{d}{\mathrm{O}}_{\mathrm{2}}}{\mathrm{d}t}}={\displaystyle \frac{F_{{\mathrm{O}}_{\mathrm{2}}\left(\mathrm{s}\right)}-F_{{\mathrm{O}}_{\mathrm{2}}\left(\mathrm{e}\right)}}{h}}+{\mathrm{adv}}_{{\mathrm{O}}_{\mathrm{2}}},\end{array}$$

where $F_{{\mathrm{O}}_{\mathrm{2}}\left(\mathrm{s}\right)}$ is the net surface O₂ flux at the canopy, $F_{{\mathrm{O}}_{\mathrm{2}}\left(\mathrm{e}\right)}$ is the O₂ entrainment flux, h is the boundary layer height, and adv${}_{{\mathrm{O}}_{\mathrm{2}}}$ is the advection term. The surface flux is calculated with Eq. (10) and the entrainment flux is based on the following equation (see also Eq. 2):

$$\begin{array}{}\text{(A2)}& F_{{\mathrm{O}}_{\mathrm{2}}\left(\mathrm{e}\right)}=-w_{\mathrm{e}}\cdot {\mathrm{\Delta }}_{\text{(ft-bl)}}{\mathrm{O}}_{\mathrm{2}},\end{array}$$

where w_e is the entrainment velocity and Δ_(ft-bl)O₂ is the jump of O₂. The jump of O₂ was determined the same way as for CO₂, by tuning the initialization of the jump until the decrease or increase of CO₂ $/$ O₂ of the model matched with the observational data.

A3 Validation of CLASS

Figures A3 and A2 present a comparison between the model output of CLASS and the corresponding measurements for the representative days of 2018 and 2019, assessing various parameters. Both figures demonstrate that the model compares well to the observed data. CLASS accurately follows the observed temperature increase (Fig. A2a). A constant difference of approximately 8K between 2018 and 2019 is seen for both the model and the observations. This persistent difference is attributed to a heat wave rather than a drought in Hyytiälä, as a drought would have intensified the divergence between the 2018 and 2019 simulations throughout the day. Moreover, CLASS adequately models specific humidity for both years, assuming an initial relative humidity of 80 % for 2018 (Fig. A2b). The sensible heat flux (Fig. A2c) and latent heat flux (Fig. A2d) exhibit minimal differences between the 2018 and 2019 simulations. The accurate representation of atmospheric properties in CLASS consequently results in a satisfactory representation of boundary layer height development for both years in comparison to the observed data from radiosondes (Fig. A3a)

The various CO₂ fluxes simulated by CLASS exhibit a high level of agreement with the observational data for both 2018 and 2019 (Fig. A3d and e). While there are subtle differences evident between the observations for the two years, CLASS adeptly captures these nuances. Consequently, the model provides an accurate representation of plant behaviour under both normal and warmer conditions. The elevated temperatures (+8 K) and slightly reduced soil moisture (−0.03 m³ m⁻³) contribute to a slightly higher GPP and TER flux. Our study reaffirms that the vegetation in Hyytiälä did not undergo any stress during the 2018 European drought, which would have resulted in a lower GPP and lower latent heat flux (Lindroth et al., 2020).

For the 2018 case, we altered only a few initial conditions (see Table C3). However, both the decrease in CO₂ and the increase in O₂ during the day exhibit close similarity between the model and the observations. This outcome underscores that even with minimal changes in the initial conditions for the 2018 case and while keeping the other variables constant (e.g. the jumps), we can successfully replicate a realistic new day based on the base case.

It is important to note that only the net ecosystem exchange (NEE) data are obtained directly from eddy covariance measurements. The gross primary production (GPP) is inferred from a light- and temperature-based function and the total ecosystem respiration is calculated as the residual between NEE and GPP (Kulmala et al., 2019; Kohonen et al., 2022). This distinction may explain the challenge in aligning the TER flux of the observations with the model, as the model exhibits notable discrepancies from the observations for the 2018 and 2019 cases. The model's simulated respiration increase based on temperature appears more extreme compared to the observations. However, several studies (Lindroth et al., 2008; Gao et al., 2017; Heiskanen et al., 2023) indicate that the model's increase in TER between 2018 and 2019 is slightly too high, while the change based on observations is too low. As a result, it is plausible that the true respiration flux lies somewhere between the model output and the observational data.

Appendix B: Figures

Figure B1 Similar to Fig. 5 but now for the base case (2019) and the background sensitivity studies with a lower jump ratio between O₂ and CO₂ (lower Δ_(ft-bl)), and with a lower initial jump for CO₂ (lower initial Δ_(ft-bl)). The diurnal variability of the exchange ratio of the atmosphere is now added (ER_atmos: (h)).

Figure B2 Similar to Fig. 5 but now for the base case (2019) and the dry and warm sensitivity studies with a high soil moisture and a low soil moisture, both with higher air temperatures compared to the 2019 base case. The diurnal variability of the exchange ratio of the atmosphere is now added (ER_atmos: (h)).

Appendix C: Tables

Table C1 The initial conditions used for the three sensitivity analyses compared to the initial conditions for the 2019 base case. The subscript (0) indicates the first time step.

Table C2 Initialization of the CLASS model for the 2019 base case based on 10 July 2019. The initialization is based on the SMEAR II data (Hari et al., 2013), our OXHYYGEN campaign data (radiosondes or O₂ and CO₂ measurements) (Faassen et al., 2023), and studies that show ranges for parameters for the plants and soil (Lindroth et al., 2008; ECMWF IV, 2014; Vilà-Guerau de Arellano et al., 2015).

Table C3 Adjustments for the 2018 case (warm case) compared to the 2019 values shown in Table C2. Only the initial potential temperature (θ₀), initial soil moisture (w_g), and CO₂ mole fraction (CO_2,0) are adjusted based on the aggregate of 28 and 29 July 2018. It was assumed that the initial relative humidity stayed constant at 80 % with increasing temperatures; therefore, the initial specific humidity was also adjusted.