CloudRoots: integration of advanced instrumental techniques and process modelling of sub-hourly and sub-kilometre land–atmosphere interactions

The CloudRoots field experiment was designed to obtain a comprehensive observational dataset that includes soil, plant, and atmospheric variables to investigate the interaction between a heterogeneous land surface and its overlying atmospheric boundary layer at the sub-hourly and sub-kilometre scale. Our findings demonstrate the need to include measurements at leaf level to better understand the relations between stomatal aperture and evapotranspiration (ET) during the growing season at the diurnal scale. Based on these observations, we obtain accurate parameters for the mechanistic representation of photosynthesis and stomatal aperture. Once the new parameters are implemented, the model reproduces the stomatal leaf conductance and the leaf-level photosynthesis satisfactorily. At the canopy scale, we find a consistent diurnal pattern on the contributions of plant transpiration and soil evaporation using different measurement techniques. From highly resolved vertical profile measurements of carbon dioxide (CO₂) and other state variables, we infer a profile of the CO₂ assimilation in the canopy with non-linear variations with height. Observations taken with a laser scintillometer allow us to quantify the non-steadiness of the surface turbulent fluxes during the rapid changes driven by perturbation of photosynthetically active radiation by cloud flecks. More specifically, we find 2 min delays between the cloud radiation perturbation and ET. To study the relevance of advection and surface heterogeneity for the land–atmosphere interaction, we employ a coupled surface–atmospheric conceptual model that integrates the surface and upper-air observations made at different scales from leaf to the landscape. At the landscape scale, we calculate a composite sensible heat flux by weighting measured fluxes with two different land use categories, which is consistent with the diurnal evolution of the boundary layer depth. Using sun-induced fluorescence measurements, we also quantify the spatial variability of ET and find large variations at the sub-kilometre scale around the CloudRoots site. Our study shows that throughout the entire growing season, the wide variations in stomatal opening and photosynthesis lead to large diurnal variations of plant transpiration at the leaf, plant, canopy, and landscape scales. Integrating different advanced instrumental techniques with modelling also enables us to determine variations of ET that depend on the scale where the measurement were taken and on the plant growing stage.

1 Introduction

Evapotranspiration (ET), the net exchange of water vapour between the land and the atmosphere, remains an elusive process to be measured, quantified, and represented in models because it depends on the interaction of multiple processes that act in a wide range of scales (Katul et al., 2012). ET is a key variable in the exchange of heat, moisture, and carbon dioxide at the surface, and it strongly depends on how radiation and energy are partitioned into latent and sensible heat (Moene and Dam, 2014; Monson and Baldocchi, 2014). The amounts of direct and diffuse radiation reaching the leaves depend on the transfer of radiation that is strongly perturbed by clouds and aerosols and on its subsequent penetration into the canopy. Triggered by ambient light conditions, the stomatal responses coupled to the surface and boundary layer dynamics is the main driver that regulates how the net available radiative energy is partitioned between the turbulent sensible and latent heat fluxes (van Heerwaarden and Teuling, 2014). However, due to the highly non-stationary nature of atmospheric radiation (van Kesteren et al., 2013b) and turbulent nature of the meteorological fluctuations, we still lack a fundamental understanding of the two-way feedback between stomatal control and cloud radiation perturbations across scales and land and atmosphere conditions (Katul et al., 2012; Sikma et al., 2018).

The bidirectional link between surface processes and boundary layer clouds as described above is what we refer to as the CloudRoots concept, where boundary layer dynamics and clouds are rooted in or coupled to the surface and vice versa (Vilà-Guerau de Arellano et al., 2014). The degree of coupling depends on soil, plant, and weather conditions characterized by the diurnal variability of wind, temperature, and specific humidity (Sikma et al., 2018). To fully comprehend this system requires inclusion of all necessary parameters at the required spatial scales, from the size of the stomata (10–100 µm) to the depth of the boundary layer and cloud top (∼3 km), as well as temporal scales from seconds to daily and seasonal cycles and across disciplines, bringing together experts from diverse fields from ecophysiology to turbulence. This can only be obtained by integrating experimental and modelling efforts. Here we describe and show the first results of the CloudRoots field experiment aimed at obtaining new understanding about the interaction between the soil, vegetation, and the clear–cloudy boundary layers at these sub-hourly and sub-kilometre scales, i.e. on spatio-temporal scales smaller than the characteristic grid resolution scales of the weather (typical resolution ranging from 1 to 10 km) and climate (typical resolution ranging from 20 to 100 km) models. In that respect, the CloudRoots field campaign continues the tradition of experiments that connect land surface properties with boundary layer dynamics but now instead by using advanced instrumental techniques and by modelling the coupling between the essential processes. Two examples of such previous campaigns are the First ISLSCP Field Experiment (FIFE) (Hall et al., 1989) and the Boreal Ecosystem-Atmosphere Study (BOREAS) (Sellers et al., 1995).

Thanks to their high-quality routine measurement programme (Franz et al., 2018; Rebmann et al., 2018), ICOS sites lend themselves as anchors for additional experiments. Here, we describe the CloudRoots campaign near the agricultural site “Selhausen” (ICOS site DE-RuS) and the Jülich Observatory for Cloud Evolution – Core Facility (JOYCE, http://joyce.cloud, last access: 21 August 2020) in Germany during spring 2018 (Löhnert et al., 2015). In order to quantify all the necessary scales of interest (leaf, canopy, and landscape), we complemented the existing radiation, flux, and soil measurements of the ICOS site by scintillometry, microlysimeters, sap flow and leaf-level flux measurements, quasi-instantaneous vertical profiles, and spectroscopic measurements of vegetation indices and sun-induced fluorescence (SIF). Scintillometers provided minute-scale turbulent fluxes enabling us to connect stomatal responses to the energy, moisture and carbon dioxide (CO₂) fluxes at this timescale. Microlysimeters, soil flux chambers, sap flow, leaf-level chambers, and canopy-resolving profiles all have the ability to distinguish vegetation from soil CO₂ and water vapour (H₂O) fluxes in contrast to the eddy covariance technique that provided net fluxes from the two sources combined. The remote sensing measurements of boundary layer dynamic evolution and cloud properties made at JOYCE provided evidence on diurnal variations of the boundary layer depth, the role of entrainment, and cloud diurnal variability. A key aspect of the research strategy of CloudRoots is the integration of all these measurements in a land–atmosphere conceptual model CLASS (Vilà-Guerau de Arellano et al., 2015). This model has been specially developed to support the interpretation of measurements at the sub-hourly scales (Vilà-Guerau de Arellano et al., 2019).

To this end, we study the following five facets of the diurnal interactions between the land and the atmosphere: (i) observational validation at leaf level of the mechanistic model representation of the stomatal aperture and photosynthesis, (ii) the diurnal variability of H₂O–CO₂ flux partition due to the soil and plant contributions at the canopy level, (iii) the non-steadiness of these fluxes due to the influence of clouds, (iv) the spatial heterogeneity of ET inferred from the SIF measurements, and (v) the integration of the observations in the conceptual model CLASS to quantify the influence of land–surface heterogeneity and advection. We finally obtain a daily estimation of ET and discuss differences with respect to the observational or modelling techniques.

The paper is organized as follows. In Sect. 2 we give a detailed overview of the field experiment with special emphasis on the instrumentation used that serve the overall goals of our CloudRoots concept. The results in Sect. 3 are organized into the five topics outlined below. First, at leaf level, we validate a photosynthesis–conductance mechanistic model that is commonly used in large-eddy simulations (Pedruzo-Bagazgoitia et al., 2017; Sikma et al., 2018) and the global numerical model prediction system ECMWF-IFS (Boussetta et al., 2013). This allows us to assess the need to revisit currently used constants in the mechanistic model representing photosynthesis. This part is completed by comparing leaf transpiration rate with tiller-level measurements of sap flow at different stages of the growing season. Second, and in order to scale up to the canopy level, we analyse the soil and plant partitioning of the net ET and net ecosystem exchange (NEE) based on the inversion of observed high-resolution vertical concentration profiles (Warland and Thurtell, 2000; Santos et al., 2011). Third, in analysing the impact of clouds on ET, we measure the potential effectiveness of diffuse radiation in enhancing ET and NEE (Kanniah et al., 2012). Extending previous work by van Kesteren et al. (2013b), we quantify the time lag between fluctuations in incoming shortwave radiation and ET in the field. These real-world measurements are an essential addition to time lag of plant responses to radiation changes studied in laboratory experiments (Vico et al., 2011). Fourth, we infer the spatial variability of ET around the CloudRoots site using SIF remote sensing observations. Fifth, all of these observations are then integrated into several numerical experiments made by CLASS with special emphasis on the treatment and role of how to include surface heterogeneity and heat and moisture advection to improve the interpretation of the observations. Finally, in the discussion in Sect. 4 we bring together and discuss all CloudRoots methodologies by comparing their daily ET estimates. Conclusions are given in Sect. 5.

2 Description of the CloudRoots field experiment and modelling effort

2.1 Site description

The CloudRoots field campaign was carried out at the Terrestrial Environmental Observatory (TERENO) Selhausen, which is located in the southern part of the lower Rhine embayment in western Germany (50^∘52^′09^′′ N, 6^∘27^′01^′′ E, 104.5 m altitude) in a region largely dominated by agriculture (Fig. 1). In 2011, the site was equipped with micrometeorological measurement devices for long-term monitoring of energy and carbon exchange. Since 2015, the station has been extended in accordance with ICOS standards for Level 1 sites (ICOS site code DE-RuS) (Ney et al., 2020). For this campaign, a further IRGASON eddy covariance (EC) system with an open-path gas analyser (see Sect. 2.3.7) was placed on the test field and used for additional flux measurements presented here.

Figure 1 (a) Aerial view (Bing Maps, © 2019 Microsoft Corporation © 2019 DigitalGlobe © CNES (2019) Distribution Airbus DS) of the observation area. The ICOS Selhausen test site is located in the middle of the 10 × 10 km map section. The surrounding agricultural area was classified into the categories bare soil (including “late crops” following Table 3) and vegetated (“early crops”, forest and grassland following Table 3) during the IOP 1. (b) Corresponding sensible heat flux (H) during IOP 1, whereby H of bare soil and vegetated area were measured and the regional average was estimated as weighted average (60 % and 40 % for vegetated and bare soil, respectively). (c) Schematic sketch of horizontal (red) and vertical (black) length scales influencing the measurements. The larger indicated horizontal and vertical scales indicate the spatial scales of boundary layer dynamics. Horizontally, the 100 m scale is the size of the field hosting the ICOS test site.

The test field covered 9.8 ha and was surrounded by other croplands (Ney and Graf, 2018). As Fig. 1 shows, these cultivated areas are mainly comprised of winter wheat, winter barley, sugar beet, rapeseed, maize, potatoes, and peas, whereby the various field sizes and locations of crops has led to small-scale heterogeneity in the vegetation cover. An agricultural road, mainly used by farm machinery, passes by the northern edge of the field. The next inhabited settlement is located 500 m to the west (Fig. 1a). There are two lignite open-cast mines in the wider surrounding of the study site, located 6 km northeast (extension of 4400 ha with a maximum depth of 470 m b.g.l.) and 6 km west (extension of 1400 ha with a maximum depth of 200 m b.g.l.). In general, the land surface at the study site is flat and has a slope less than 4^∘. A loess layer with a thickness of about 1 m covers Quaternary sediments, which were mainly built up from fluvial deposits of the Rur river system. The overlying soil is an Orthic Luvisol according to the USDA classification (IUSS Working Group WRB, 2006), whose texture is silt loam with a mixture of 20 % clay, 67 % silt, and 13 % sand.

The local climate is classified as temperate maritime with an annual mean air temperature of 10.3 ^∘C and an annual mean precipitation of 718 mm (reference period 1981–2010, with data taken from the DWD climate station of the Forschungszentrum Jülich 5.3 km from the test site). The observation period from the beginning of May until the end of June 2018 was characterized by a 2.9 ^∘C higher mean air temperature (17.5 ^∘C) and 46 % less precipitation in comparison to the long-term average. Figure 1b shows the heterogeneity quantified by the sensible heat fluxes measured at the CloudRoots site and a bare soil field nearby. In consequence and as shown by Fig. 1c, in CloudRoots we aim to integrate horizontal and vertical scales in the analysis of ET and its relation to boundary layer dynamics.

The field campaign covered the main growing phases (booting, heading, and maturity stages) of winter wheat. During the observation period, we did three intensive observation periods (IOP). During these IOPs the following complementary instruments and measurements were added: microlysimeters, leaf-level measurements, SIF measurements on canopy and regional scale, and vertical profiles of state variables and CO₂ within and above the canopy were performed. Figure 2 shows a timeline of the deployment of the campaign-specific measurement setup (see Sect. 2.2 and 2.3) that includes the IOPs on 7 May (IOP 1), 15 June (IOP 2), and 28 June 2018 (IOP 3). The main meteorological and biometric conditions are summarized in Table 1. The test field was cultivated with a crop rotation cycle typical of the region (Ney et al., 2020). The rotation prior to the observation period was beet, potatoes, and winter wheat (catch-crop) and sugar beet. Residues of the harvest of sugar beet were left on the site and ploughed in before the cultivation cycle started with the sowing of winter wheat (Triticum aestivum L.; variety Premio) in October 2017. The field was fertilized with mineral nitrogen (N) once in March, April, and May 2018 (81.6, 39.2, and 50 kg N ha⁻¹, respectively). The wheat was harvested on 17 July 2018 with a yield of 92 dt ha⁻¹. A detailed overview of the field management practices before, during, and after the campaign is given in the Appendix (Table A1).

Figure 2 Campaign-specific measurement setup and temporal developments from May to June 2018, including three intensive operation periods (IOP).

Table 1 Meteorological and biometric conditions during the intensive operation periods on 7 May (IOP 1), 15 June (IOP 2), and 28 June 2018 (IOP 3). Global radiation, water vapour–pressure deficit (VPD), photosynthetically active radiation (PAR), and soil water content (SWC) are daily averages. The meteorological variables were measured at the height 2.4±0.1 m (see Sect. 2.3.7 for details).

^* Daily averages calculated from sunrise to sunset.

2.2 Weather and crop description during the IOPs

The weather situation during all three IOPs was mainly characterized by an anticyclonic pressure pattern over central Europe (IOP 1 and IOP 2), extending up to northern Europe during IOP 3, which led to high 2 m temperatures up to 24 to 26 ^∘C during IOP 1 and IOP 2 and 28 ^∘C during IOP 3 (Table 1). Cloudiness and temperature inversion heights at the top of the atmospheric boundary layer were different. While weak subsidence motions during IOP 1 led to a slightly rising temperature inversion layer between 1200 and 2000 m above ground level (a.g.l.) with clear conditions during the whole period (mean daytime global radiation S↓ of 514 W m⁻²), a weak cold front passed the measuring site from the northwest in the early morning of IOP 2 (mean daytime S↓ of 311 W m⁻²). Diurnal heating caused the replacement of a layer of stratocumulus at a height of 1800 m a.g.l. in the morning, followed by the appearance of scattered towering cumulus clouds. Light showers occurred only in the vicinity of the site. During IOP 3, a few shallow cumulus and cirrus clouds appeared, despite the existence of a small upper-air low that passed the area around the edge of a larger cut-off, although it was located above southeastern Europe. The mixed boundary layer was topped at a height of around 1700 m a.g.l.

The persistent high-pressure weather conditions resulted in a drought during the entire observation period. Ongoing dryness led to a reduction in the soil water content at 20 cm depth (Table 1) from 27 vol % during IOP 1 to 15 vol % at IOP 3. Maturity occurred 14 d earlier than in previous years. The leaf area index (LAI) ranged from 4.5 m² m⁻² (green growing stage) in IOP 1 to 5.5 m²m⁻² in IOP 2 (green/yellow ripening stage). No changes in LAI were observed between IOP 2 and IOP 3 (yellow senescence stage).

Table 2 List of symbols, description, units, and the representative scale.

2.3 Instrument description

Table 2 summarizes all the variables measured and modelled during CloudRoots, together with specific nomenclature and information on units and scales.

2.3.1 Microlysimeters

For direct measurements of soil evaporation (E_lys), four microlysimeters were installed at a number of locations around the EC-station (one in each cardinal direction) at the beginning of every observation period. In order to obtain an undisturbed soil monolith for each microlysimeter, an SDR-35 polyvinyl chloride (PVC) collar with an inner diameter of 0.2 m, a wall thickness of 0.005 m, and a depth of 0.11 m was pushed carefully into the ground. Afterwards the collar including the soil column was retrieved, its outside was cleaned, and the bottom of each lysimeter was sealed with an acrylic glass disc, which prevented percolation and capillary rise from or into the microlysimeter. The microlysimeters were then weighed initially and returned to their original positions. We made sure that the lysimeters were levelled with the soil surface, their walls were fully surrounded by soil, and that the crop was affected and destroyed as little as possible, so that the general conditions and characteristics of the field site could still be maintained (e.g., regarding heat flux, shading). All four microlysimeters were subsequently collected, cleaned, weighed, and distributed again every 60 or 90 min. A scale with a precision of 0.1 g (equivalent to 0.00318 mm evaporation) was used. The scale was enclosed in a box to avoid wind effects during the measurements. Finally, the measured weight differences were converted to the standard units of the surface energy balance (W m⁻²) by means of the lysimeters surface area, the time periods between weighing, and the latent heat of vaporization (Quade et al., 2019).

2.3.2 Soil CO₂ flux chambers

Soil respiration (R_s) was observed with an automated soil CO₂ gas flux system (Li-8100, Li-Cor Inc. Biosciences, Lincoln, Nebraska, USA), connected to four long-term soil flux chambers. The chambers were installed close to the EC-station (one in each cardinal direction) on top of PVC soil collars with a diameter of 0.2 m and a total height of 0.07 m, from which 0.05 m was inserted into the soil. Each chamber was closed at 30 min intervals for 90 s during flux measurements, while CO₂, water vapour concentration, and chamber headspace temperature were recorded at a sampling rate of 1 Hz. The CO₂ concentration was standardized to dry air and a constant temperature to eliminate effects of changes in air density and water vapour dilution during closure time. R_s was subsequently calculated by adjusting a linear regression fit to the final 60 s of the measurement before reopening.

2.3.3 Leaf-level measurements

Leaf gas exchange was measured using a Li-Cor LI-6400XT portable photosynthesis system with a 6400-02B LED light source. Leaf-level measurements included instantaneous stomatal conductance to water vapour (g_sw) and photosynthesis (A_leaf), maximum light-saturated photosynthesis (A_max) is the maximal primary productivity under high light conditions, CO₂-response curves, and light-response curves. Measurements of g_sw and A_leaf were performed during the three IOPs, starting at sunrise and ending when measurements of g_sw indicated that stomata had nearly closed (g_sw<0.05 mol m⁻² s⁻¹). For measurements of g_sw and A_leaf, tillers were picked randomly in the field and immediately mounted in the leaf chamber for measurements. Initial tests showed no difference in g_sw between excised and attached tillers. Settings of leaf chamber photosynthetically active radiation (PAR) and CO₂ followed the diurnal variability measured in the field. For comparison with other observations, measurements of g_sw and A_leaf were binned and averaged at 30 min intervals. A_max was measured during the three IOPs, as well as on 8 May between 10:00 and 12:00 UTC. For measurements of A_max, the light intensity (PAR) was set to 1500 µmol m⁻² s⁻¹ and the leaf was equilibrated under a reference CO₂ concentration of 450 µmol CO₂ mol⁻¹ air. CO₂ response curves were measured during IOP 1 and IOP 3, prescribing CO₂ concentrations in the following order: 0, 50, 100, 150, 250, 350, 450, 600, 800, and 1200 µmolCO₂ mol⁻¹ air. All CO₂-response curves were measured using a light intensity (PAR) of 1500 µmol m⁻² s⁻¹. Light-response curves were measured on IOP 1 only and used a reference CO₂ concentration of 450 µmol CO₂ mol⁻¹ air. PAR values were changed in the following order: 0, 25, 50, 100, 200, 400, 800, 1200, and 1500 µmol m⁻² s⁻¹. The stomatal conductance to water vapour (g_sw mol m⁻² s⁻¹) of the curves that relate A and PAR (A–PAR) in between 0 and 200 µmol m⁻² s⁻¹ for the three repetitive experiments within the PAR range were (average and standard deviation in brackets): 0.49 (0.13), 0.12 (0.02), and 0.34 (0.06). Leaves were allowed to equilibrate to leaf chamber conditions in terms of gas exchange (approximately 1–2 min) but not in terms of stomatal aperture. For all measurements, leaf chamber temperature was set between 20 and 25 ^∘C. Relative humidity in the leaf chamber was set between 60 % and 75 %. Measurements of A_max, CO₂-response curves and light-response curves were performed on attached tillers.

2.3.4 Sap flow

Sap flow in wheat tillers was measured with the heat balance method (Sakuratani 1981; Baker and van Bavel, 1987). A total of 24 tillers were selected at random, diameters measured with an electronic calliper, and SGA3-type sap-flow sensors installed at the lowest possible internodes following the procedure recommended by the manufacturer (Dynamax, 2007, 2017). Sensors were connected with electrically shielded wire to AM 16/32 multiplexers controlled and scanned by CR1000 data loggers (Campbell Scientific, Logan, Utah, USA). Energy supply to the stem heaters was carefully regulated to the highest permissible level in order to obtain a strong heat signal. We employed the dual voltage regulators (Dynamax AVRDC) which were parts of wired measurement, control, and extension units assembled and tested by the heat balance sensor manufacturer (Flow32 1K A and B models, Dynamax Inc., Houston, Texas, USA). Data were processed according to the calculation procedure of Dynamax (2007) with adaptations to wheat (Langensiepen et al., 2014) to obtain reliable data on the convective stem heat flow generated by sap flow. Here we take the evolution of the tiller densities from 480 tillers per square metre (IOP 1 and IOP 2) to 370 tillers per square metre (IOP 3) into account.

2.3.5 Profiling elevator

Vertical profiles H₂O and CO₂ expressed as mole fractions χH₂O and χCO₂ (mole of substance per mole of moist air), temperature (T_air,p) and wind speed (u_p) from the soil surface to the surface layer above the crop canopy were measured with a portable elevator system. The elevator moved continuously up and down the measuring sensors attached to an extension arm over a total profile height of 2 m. A sampling tube connected to a differential gas analyser (LI-7000, Li-Cor Inc. Biosciences, Lincoln, Nebraska, USA) collected χH₂O and χCO₂ at a frequency of 20 Hz. T_air,p and u_p were measured at the same frequency by a ventilated fine-wire thermocouple (FW3, Campbell Scientific, Logan, Utah, USA) and a hot-wire anemometer (8455-075-1, TSI, Shoreview, Minnesota, USA). All measurements were duplicated as a continuous fixed-height measurement at the top of the profile. During the data post-processing, the temporal and vertical resolution of the mean profiles was set to a time-averaging block of 30 min with a vertical resolution of 0.025 m. Time delays in each variable with respect to the position caused by response times of the sensors, electronic delays, and the tube transport of the gas samples were adjusted by a hysteresis minimization algorithm. Detailed information on the profile measurement setup and the processing of the data profile is given in Ney and Graf (2018). The measured concentration profiles were then used to determine the vertical source profiles of H₂O and CO₂, with the aim of providing an independent, non-invasive partitioning between above-ground net primary production (NPP) and R_s or evaporation (E) and transpiration (Tr_p). To estimate source profiles and flux partitioning, we used an analytical dispersion Lagrangian technique introduced by Warland and Thurtell (2000) and further developed by Santos et al. (2011). Other than in the above-mentioned literature, a simple optimization method (Nelder and Mead, 1965) was used to fit four parameters: soil source, canopy source, and shape parameters p and q of a beta distribution, which describes the vertical source distribution within the canopy.

2.3.6 Scintillometer

The receiver of a displaced-beam laser scintillometer, hereafter referred to as DBLS (SLS-20, Scintec, Rottenburg, Germany), was placed 9 m southeast of the EC station (Fig. 1). The scintillometer measurements height was 1.95 m a.g.l. The path length towards the instrument transmitter was 86.8 m. It was pointed along a transect from northwest to southeast. The DBLS measures the scintillation intensity of two displaced laser beams (wavelength of 670 nm and separation distance of ∼2.7 mm). The structure parameter of temperature ($C_T^{\mathrm{2}}$) and dissipation rate of turbulent kinetic energy (ε) are determined from the log variance of one beam and log covariance between the beams. The general equation that links the scintillometer measurements to fluxes is given by

$$\begin{array}{}\text{(1)}& F_x=\mathit{\rho }{K}_x\left(u_{*},{\displaystyle \frac{z}{L}}\right)z^{\frac{\mathrm{1}}{\mathrm{3}}}\sqrt{C_{x^{\mathrm{2}}}^{\mathrm{2}}},\end{array}$$

where F_x is defined as the turbulent flux of the transported variable x, $C_x^{\mathrm{2}}$ is the structure function parameter of x, and K_x represents the turbulent exchange coefficient that links F_x to $C_x^{\mathrm{2}}$. K_x is a function of the friction velocity, u_*, and the Obukhov length, L. Finally, ρ is the air density, and z is the measurement height above the surface. For the sensible heat flux H=F_T, x represents temperature (T) and appropriate constants need to be added to convert Eq. (1) to energy fluxes. H, u_*, and L are solved iteratively as a function of the DBLS-measured $C_T^{\mathrm{2}}$ and ε (Thiermann, 1992; Hartogensis et al., 2002). The Monin–Obukhov Similarity Theory (MOST) functions that define K_x were taken from Kooijmans and Hartogensis (2015). For our purpose, however, the exact shape of the MOST functions is of minor importance as we are primarily interested in the dynamic temporal behaviour of the fluxes rather than an accurate description of their quantitative values. We are aware that advective contributions can lead to the violation of MOST. However, advection did not influence our measurements for two reasons. First, the scintillometer transmitter and receiver are far enough from the edges of the CloudRoots field given the height of the sensor (1.95 m), the wind speed and direction during the IOPs, and the stability conditions. All of these mean that the footprints are small enough to fit within the field. The typical footprint length (90 % footprint contribution) for the 3 IOPs yields IOP 1 (85 m), IOP 2 (30 m), and IOP 3 (75 m). Second, the scintillometer has a path weighting function that is at its maximum in the middle of the path and near-zero at the transmitter and receiver positions, i.e. the major contribution occurs at the farthest point of the field edge.

The added value of DBLS fluxes over the traditional EC method is that they converge to statistically stable flux estimates at much shorter flux averaging times of 1 min or less, while the EC technique typically requires flux averaging times of 10 to 30 min (Hartogensis et al 2002; van Kesteren et al., 2013b). The essence behind this is that the flux estimate is based on structure parameters that are defined in the inertial range of the turbulent spectrum. As such, the flux estimates rely on a limited range of the turbulent scales that contribute to the flux rather than all of them as is the case with the EC method.

We also adopted the combination technique introduced by van Kesteren et al. (2013a, b) to obtain fluxes of H₂O and CO₂ at these short timescales. This technique combines structure parameters of H₂O and CO₂ that are obtained from H₂O and CO₂ time series from an Infra-Red Gas Analyser (IRGASON system; see Sect. 2.3.7) with an exchange coefficient defined by the DBLS fluxes to finally calculate flux estimates of H₂O and CO₂. In other words, with u_* and L solved with the DBLS, Eq. (1) can be evaluated using structure parameters of trace gases x, where in this case x represents the specific density, q_x, of H₂O or CO₂.

2.3.7 Eddy covariance and ancillary micrometeorological measurements

A continuously running EC system was operated in the middle of the field (Fig. 1), comprising a three-dimensional sonic anemometer (Model CSAT-3, Campbell Scientific, Inc., Logan, Utah, USA) and an open-path infrared gas analyser (Model LI-7500, Li-Cor, Inc., Biosciences, Lincoln, Nebraska, USA). The sensors height was 2.34 m a.g.l. Raw data were sampled in 20 Hz mode and fluxes and averages were calculated as 30 min block averages using the TK3.11 software package developed at the University of Bayreuth, including corrections and quality control as given in Mauder and Foken (2011). Missing values in the calculated turbulent fluxes were filled with the marginal distribution sampling (MDS) method following Reichstein et al. (2005) which is implemented in the REddyProc software package (Wutzler et al., 2018). The station also included measurements of all components of the radiation budget (NR01, Hukseflux, Delft, the Netherlands), PAR (LI-190R, Li-Cor Inc. Biosciences, Lincoln, Nebraska, USA, and BF5, Delta-T Devices, Cambridge, UK), air temperature (T_air), and humidity (HMP45C, Vaisala Inc., Helsinki, Finland) at 2.4 m, and precipitation (Thies Clima type tipping bucket, distributed by Ecotech, Bonn, Germany) at 1.0 m a.g.l. Radiation measurements were taken at 2.5 m. Soil heat flux, temperature, and moisture were measured next to the station, more specifically, this was performed using the following parameters: 3 times HFP01SC at 3 and 8 cm, Hukseflux, the Netherlands; 3 times TCAV, Campbell Scientific, Logan, USA; and 1 cm, 5 cm, and 2 to 65 cm layer average, 2 times CS616, Campbell Scientific, Logan, USA, 2 to 6 cm layer average. In addition, we measured at five points distributed across the field using the wireless SoilNet sensor system (Bogena et al., 2010). One SoilNet point was placed next to the station, while the other four were placed next to the soil CO₂ efflux chambers described above. Each SoilNet point comprised a single soil heat flux measurement at 5 cm (HFP01SC, see above) and combined temperature and soil water content measurements in depths of 1, 5, 10, 20, 50, and 100 cm (SMT100, Truebner GmbH, Neustadt, Germany).

A second mobile EC station with instruments heights of 1.93 m a.g.l. was deployed in the immediate vicinity of the continuously monitoring station during the measurement campaign. The system comprised an IRGASON EC system (SN1185 Irgason EC150, Campbell Scientific, Inc., Logan, Utah, USA; PTB101B pressure sensor, Vaisala Inc., Helsinki, Finland) with an additional LI-7500 sensor (same manufacturer). Here, fluxes were processed with the LiCor EddyPro v6.2.2 software. Radiation (CM11 for global and CG2 for long-wave radiation, Kipp & Zonen B.V., Delft, Netherlands), ground heat flux (4× HFP01SC at 5 cm depth, Hukseflux, the Netherlands), and temperatures at depths of 2 cm (4×) and 8 cm (2×) were also measured at this station.

2.3.8 Canopy-level measurements of reflectance and sun-induced fluorescence (SIF): FloxBox

A field spectroscopy system was used (FLOX, JB Hyperspectral Devices UG, Düsseldorf, Germany) for canopy-level measurements of reflectance and SIF. FLOX is constructed for high temporal frequency acquisition of continuous top-of-canopy optical properties with a focus on sun-induced chlorophyll fluorescence. The system is equipped with two spectrometers: an Ocean Optics FLAME S, covering the full range of Visible and Near-Infrared (VIS-NIR) and an Ocean Optics QEPro, with a high spectral resolution (Full Width at Half Maximum, FWHM, of 0.3 nm) in the 650–800 nm range of the fluorescence emission. The optical input of each spectrometer is split between two fibre-optic cables that lead to a cosine receptor that measures solar irradiance and a bare fibre bundle that measures the target-reflected radiance. Spectrometers are housed in a Peltier thermally regulated box to keep the internal temperature lower than 25 ^∘C in order to reduce dark current drift. The signal is automatically optimized for each channel at the beginning of each measurement cycle, and two associated dark spectra are collected as well. Metadata such as spectrometer temperature, detector temperature and humidity, Global Positioning System (GPS) coordinates, and time are also simultaneously stored in the secure digital memory of the system. More detailed information about the system can be found in Wohlfahrt (2018) and in Campbell (2019).

2.3.9 Regional level measurements of reflectance and sun-induced fluorescence (SIF): HyPlant

An airborne high-performance imaging spectrometer (HyPlant) was used for regional-level measurements of the same quantities. Several flight lines over the 15 km × 15 km study site with 1–3 m pixel resolution were used. HyPlant is a hyperspectral imaging system for airborne and ground-based use, developed as a cooperative effort between Forschungszentrum Jülich GmbH (Germany) and the company SPECIM (Oulu, Finland). It consists of two sensor heads, named DUAL and FLUO. The DUAL module is a line-imaging push broom hyperspectral sensor, which provides contiguous spectral information from 370 to 2500 nm in a single device that utilizes a standard objective lens with 3 nm spectral resolution in the VIS-NIR spectral range and 10 nm spectral resolution in the SWIR spectral range. The FLUO module measures the vegetation fluorescence signal with a separate push broom sensor that produces data at high spectral resolution (0.25 nm) in the spectral window between 670 and 780 nm. The position and altitude sensor (GPS/INS sensor) provides, synchronously with the image data, aircraft position and altitude data for image rectification and geo-referencing. Both imagers are mounted in a single platform with the mechanical capability to align the field of view (FOV). A more detailed description of the sensor is given in Rascher et al. (2015).

Sun-induced fluorescence (F₆₈₇ and F₇₆₀) was retrieved in the two oxygen absorption bands according to the improved Fraunhofer Line Discrimination (iFLD) method. The iFLD method was initially proposed by Alonso et al. (2008) and was adapted to allow SIF retrievals from the FLUO module of the HyPlant sensor (Rascher et al., 2015). Surface reflectance and vegetation indices were calculated after an atmospheric correction using the MODTRAN software package was applied. The atmospheric correction was performed using the MODTRAN software package (for an overview of the data processing of HyPlant, see Siegmann et al., 2019). For the reason of easier comparison of SIF values with other methods of this paper, the commonly used SIF units (mW m⁻² sr⁻¹ nm⁻¹) were replaced by a substitute (nmol m⁻² sr⁻¹ s⁻¹) using conversion factors of 6.35 for and 5.74 for F₆₈₇ and F₇₆₀, respectively.

2.3.10 Boundary layer and cloud remote sensing measurements

The JOYCE remote sensing facility (Löhnert et al., 2015) (located at a distance of 5 km from the test site) provided continuous information about boundary layer and cloud characteristics. Specifically, microwave and lidar measurements were used to compare the CLASS model results (see Sect. 2.4) with the inferred boundary layer depth. This comparison was completed by vertical profiles measured by the routine radio soundings at Essen (station ID EDZE/10410 at a distance of 75 km).

2.4 Modelling from leaf to landscape scales: CLASS

The Chemistry Land-surface Atmosphere Soil Slab (CLASS, https://classmodel.github.io/, last access: 21 August 2020) is a model that couples the soil–vegetation–atmospheric processes and is used to interpret the observations and analyse the interaction of scales (Vilà-Guerau de Arellano, et al., 2015). It contains a leaf-level representation of photosynthesis and stomatal aperture (leaf resistance). By upscaling this leaf resistance to the canopy level (surface canopy resistance), it connects with the soil processes and boundary layer diurnal dynamics. In Sect. 2.4.1 and 2.4.2, we will subsequently discuss the two main modules of CLASS that we will target in this paper, i.e. the leaf level photosynthesis module and the mixed-layer module.

2.4.1 Modelling leaf-level photosynthesis

Leaf-level photosynthesis was modelled using the representation of photosynthetic biochemistry, as included in CLASS (Vilà-Guerau de Arellano et al., 2015), which was originally developed by Goudriaan (1986) and further adapted to meteorological applications by Jacobs and de Bruin (1997). As this model describes the relationship between stomatal conductance (g_s) and photosynthesis (A), it is usually referred to as the A-g_s sub-model. In short, plant transpiration and CO₂ assimilation as part of the surface energy balance model are represented by a two-big-leaves model, one for sunlit leaves and one for shaded leaves (Jacobs and de Bruin, 1997; Pedruzo-Bagazgoitia et al., 2017). The exchange at the leaf surface depends on the gradient of atmospheric CO₂ and an internal leaf CO₂ concentration, which depends on the water vapour deficit and leaf conductance. The CO₂ exchange is upscaled to the canopy level by integrating over the leaf area index (LAI).

Available field measurements were used for improving the model settings at the leaf level. The parameters representing the initial value of the light use efficiency (α₀) and the temperature-normalized maximum leaf-level photosynthesis rate (A_m,max298) were fitted using light-response curves (Fig. 5) and CO₂-response curves (Fig. 3b) collected on 8 May 2018 (1 d after IOP 1). Table 3 summarizes the optimized values used in the A-g_s (sub)model to simulate the leaf-level photosynthesis. The A–PAR curves contain only the lower light intensity values (0–200 µmol m⁻² s⁻¹) for which the light response is near linear and not limited by CO₂ diffusion into the leaf. As leaf-level measurements of A_max indicated a decline in photosynthetic capacity in the course of the growing season (Fig. 5c), we performed additional measurements of A_m,max298 to represent the observed seasonal decline for IOP 2 and IOP 3. The impacts on these optimized values are shown and discussed in the Supplement.

2.4.2 Modelling the diurnal variability of landscape surface fluxes and boundary layer dynamics

The fundamental assumption of the mixed-layer model is that under convective conditions the atmospheric boundary layer (ABL) dynamics lead to profiles of the meteorological state variables that are uniform (well-mixed) with height. As a result, these state variables are governed by horizontally averaged 0-dimensional slab equations: one equation for the evolution through time of the slab variable and another for the difference between the residual layer (in the morning transition) and the free tropospheric values and the slab value, i.e. the jump at the interface between residual layer and ABL. The ABL dynamics are governed by the mixed-layer equations of potential temperature (heat), specific humidity (moisture), CO₂, and two horizontal wind momentum components. In addition, there is an equation that governs the boundary layer growth, which depends on the buoyancy flux at the surface and the jump in the virtual potential temperature at the interface between the atmospheric boundary layer and the free troposphere.

A key feature of the model is its representation of the sub-daily variability of the land–atmosphere interactions (van Heerwaarden et al., 2010; Vilà-Guerau de Arellano et al., 2015). The net ecosystem exchange is calculated as a result of the assimilation of CO₂ by plants and the CO₂ soil efflux. We calculate the assimilation rate from photosynthesis and the stomatal aperture measurements at leaf level (see Sect. 2.4.1), upscaled to canopy level (Ronda et al., 2001). This model depends on the diurnal variability of PAR, temperature (T_air and T_air,p), and the water vapour deficit (VPD). The two-big-leaves approach is used (sunlit and shaded) to take the different contributions of direct and diffuse radiation into account (Pedruzo-Bagazgoitia et al., 2017). The soil efflux is calculated as a function of the soil temperature and moisture. Other relevant physical processes include a radiation transfer model, the Penman–Monteith equation included in the surface energy balance, and the possibility of adding large-scale forcings such as vertical subsidence motions and large-scale advection of momentum, heat, moisture, and CO₂. Within the context of CloudRoots, it is important to mention that the model assumes a horizontal homogeneous surface. While the experimental field itself is quite homogeneous, it is surrounded by other land use types at a spatial scale that will affect the boundary layer. In that respect, and in setting the initial and boundary conditions for the numerical case, we assume that the boundary layer dynamic is governed by a sensible heat flux that is an aggregate of all the fields shown in Fig. 1b.

3 Results: integrating spatio-temporal scales from leaf to boundary layer

This section is structured following the five facets of the diurnal interactions between the land and the atmosphere outlined in the introduction.

Table 3 Parameters representing the maximum leaf-level photosynthesis rate (A_m,max298) and the initial value of light-use efficiency (α₀) under low light, as adjusted in the original A–g_s model to represent plant-specific photosynthesis characteristics for winter wheat (ww). A_m,max298 was initially fitted using the A–Ci curves and α₀ is fitted using the A–PAR curves taken during IOP 1 (Fig. 5). For IOP 2 and IOP 3, A_m,max298 values were fitted only on leaf-level measurements of A_max. The values of IOP 1 were used as numerical settings for the CLASS model runs (Fig. 16). The equivalence to typical values of the commonly used in the Farquhar–Berry–von Caemmerer (FBvC) model of leaf photosynthesis (Farquhar et al., 1980) is given in Table S1 at the Supplement.

3.1 Leaf-level exchange of H₂O and CO₂: observations and modelling

We combine leaf-level and sap flow measurements of tiller assimilation and transpiration with leaf-level assimilation modelled by CLASS, A-g_s representation, to study their variation during the growing season and the impact of unsteady PAR due to the presence of clouds.

3.1.1 Stomatal conductance and sap flow

Our leaf-level measurements revealed clear diurnal patterns in g_sw during all the IOPs (Fig. 3). The observed daily maximum g_sw decreased over the growing season. This daily maximum g_sw occurred at an earlier time during each IOP. Specifically, the 30 min average daily maximum g_sw declined from 0.84 mol m⁻² s⁻¹ (around 10:00 UTC, 12:00 LT) during IOP 1 and 0.83 mol m⁻² s⁻¹ (around 10:00 UTC) during IOP 2 to 0.30 mol m⁻² s⁻¹ (between 05:30 and 06:30 UTC) during IOP 3. The weather during IOP 2 was characterized by large cumulus clouds passing over the field site, which were made visible in the large fluctuations in PAR (Figs. 3b, 11 and 12). The cloud-related changes in light intensity induced consistent stomatal opening–closing responses during IOP 2. The relatively low g_sw observed during IOP 3 probably reflects the continuing drought that characterized the 2018 growing season in combination with the relatively high VPD and high temperatures. Sap flow measurements were performed during IOP 2 and IOP 3 (Fig. 3b and c) and one earlier non-IOP day (7 June) (Fig. 4). Measurements of sap flow revealed clear diurnal patterns for all measurement days and consistent responses to cloud-induced changes in light intensity during IOP 2 (Fig. 3b). These responses were comparable to the observed responses in g_sw during IOP 2. Interestingly, the notable decline in leaf-level g_sw between IOP 2 and IOP 3 was neither reflected in the measurements of sap flow nor in the ET measurements with the eddy covariance. For IOP 3, the ET measured by the eddy covariance had still maximum values of 300 W m⁻². Thereafter, the decrease in ET started 1 week afterwards (5 July) with values lower than 100 W m⁻². This discrepancy could partly be explained by increases in VPD and wind speed between IOP 2 and IOP 3. The more probable causes are senescence effects on physiological control of transpiration and the physical reactions to heat of the wheat tillers, which were noticeably wilting between IOP 2 and IOP 3. This observation has so far not been reported in the literature. Further studies of the relationships between senescence and simultaneously occurring changes in the heat and physical properties of wheat tillers are needed to explain this phenomenon.

Figure 3 Diurnal changes in photosynthetically active radiation (PAR) and vapour pressure deficit (VPD) measured for (a) IOP 1, (b) IOP 2, and (c) IOP 3. Leaf-level measurements of stomatal conductance of water vapour (g_sw) in (b) and (c) compared to tiller-level measurements of sap flow (E_sap). Leaf-level measurements of g_sw (blue markers) were averaged over 30 min intervals (blue line). Sap flow measurements represent the 1 SD confidence interval (shaded region) of measurements on 24 tillers averaged over 30 min timescales.

3.1.2 Observed versus modelled leaf-level photosynthesis

One of the main aims in CloudRoots is to improve the mechanistic modelling of photosynthesis and stomatal aperture. To this end, we calibrate the constants of the A-g_s model using systematic in situ field observations. Figure 5 shows the dependencies of leaf-level photosynthesis of A_leaf on PAR (Fig. 5a) and the leaf-internal CO₂ concentration (Fig. 5b) and the long-term decline in maximum light-saturated photosynthesis (Fig. 5c). Our observations indicate the need to calibrate the model depending on the functional type of the plant, in particular the dependence of A_leaf on PAR, during the field campaign. Table 2 summarizes the new constant values used in the A-g_s model adjusted to the winter wheat crop conditions.

Figure 6 shows a comparison of the model results of A_leaf using the new constants and the measurements of A_leaf and NPP together with the diurnal variation in PAR and VPD during the three IOPs. Our measurements and model results of A_leaf showed clear diurnal patterns during each IOP, and a consistent decline over the three IOPs. The decline in A_leaf was comparable to the decline in A_max (Fig. 5c) and probably reflects a combination of seasonal decay in photosynthetic capacity and increasing stomatal limitations owing to persistent drought, especially during IOP 3. The magnitude of the seasonal decline in A_leaf was comparable to the seasonal decline in NPP derived from EC data. Cloud-induced changes in PAR during IOP 2 also induced changes in A_leaf. The A-g_s model reproduced the diurnal patterns in A_leaf during each IOP, as well as the cloud-induced changes in A_leaf during IOP 2. The agreement is very satisfactory during IOP 1, which was characterized by cloudless conditions and the maturity of winter wheat. The model underestimated A_leaf during IOP 3, which was a result of the strong stomatal limitations that influenced the measurement of A_max on which the model parameterization from IOP 3 was based. The model furthermore overestimates the decline in A_leaf between 14:00 and 19:00 UTC, which probably reflects a misrepresentation of the temperature and VPD sensitivity of Triticum aestivum.

Figure 4 Sap flow measured using the heat balance method for 7 June 2018 (non-IOP day).

3.2 Canopy-level partitioning of the net H₂O and CO₂ fluxes between soil and plant processes

Moving from leaf to canopy scale, we analyse the detailed profiles of micrometeorology and carbon dioxide collected using the elevator and infer vertical assimilation profiles and the diurnal variability in the surface contributions to ET and NEE.

Figure 5 Measurements of leaf-level photosynthesis (A_leaf) as a function of photosynthetically active radiation (PAR) (a) and leaf interior CO₂ concentrations (c_i) (b). These results were used to parameterize the A-g_s model for IOP1, as indicated by the black line and shaded 1 SD confidence interval. The red line indicates the model response using the default parameter values. (c) Observed and modelled seasonal decline in maximum light-saturated photosynthesis (A_max). Boxes indicate the variability in observed values in A_max; red markers indicate the modelled net photosynthesis rate using fitted values for A_m,max298. Fitted and default A-g_s model parameter values are indicated in Table 3.

Figure 6 Measured leaf-level photosynthesis (A_leaf) compared to modelled A_leaf using the A-g_s model and canopy-level net primary productivity (NPP_canopy) for (a) IOP 1, (b) IOP 2, and (c) IOP 3. Measurements of A_leaf were plotted as 30 min averages (blue line) and their 1 SD confidence interval (shaded region). Panels (d), (e), and (f) show diurnal changes in photosynthetically active radiation (PAR) and vapour pressure deficit (VPD) measured for each IOP.

3.2.1 Concentration profiles of H₂O and CO₂, temperature and wind speed

Figure 7 shows selected 30 min mean profiles of χH₂O and χCO₂, temperature, and wind speed versus height (z) above ground level during IOP 1 and IOP 2. Over the diurnal cycle, χCO₂ concentrations fell between 08:00 and 13:00 UTC from 370 to 360  µmol mol⁻¹ in the mid-canopy during IOP 1 but stagnated slightly below 370 µmol mol⁻¹ during IOP 2. This seasonal reduction in CO₂ uptake was also observed in measured A_leaf, i.e. see the decrease of the maximum values in Fig. 6. The lowest values were observed during local noon, simultaneous with the highest PAR values (Fig. 5b). χCO₂ minima were located in the upper third of the canopy during IOP 1 and during the middle third during IOP 2. The highest χCO₂ values were found near the soil surface due to soil respiration, lower light intensity caused by shadowing, and a low amount of photosynthetic organs in the stems. Maximum χCO₂ concentrations were measured in the morning and evening hours and peaked at about 475 and 420 µmol mol⁻¹ during IOP 1 and IOP 2, respectively. The photosynthetic CO₂ uptake by plants is highly related to plant transpiration. Consequently, χH₂O in the canopy space was higher than in the air above the canopy. The highest values were found directly above the soil surface and were caused by evaporation and due to plant transpiration within the canopy.

Figure 7 Selected (08:00, 13:00 and 18:30 UTC) 30 min mean profiles of the H₂O and CO₂ mole fractions (χH₂O, χCO₂), wind speed (u_p), and temperature (T_air,p) measured at high vertical resolution during IOP 1 (upper row) and IOP 2 (lower row). Shaded areas indicate the 95 % confidence interval resulting from the standard deviation between individual profiles sampled within a 30 min average interval. The dashed lines indicate the canopy heights.

The highest temperatures appeared near the canopy top (Fig. 7d, e, j and l). In the late morning of IOP 2, the temperature reached a distinct maximum just below the canopy top (Fig. 7j). This phenomenon has been reported in previous studies (Ney and Graf, 2018) and is caused by the changing solar incidence angle. A low angle of incidence in the morning and afternoon limited the heating to an area just below the canopy surface. Previous studies have shown that the presence of such a pronounced temperature maximum has the potential to increase thermal stability within the canopy and thus inhibit the vertical turbulent exchange of sensible heat (Gryning et al., 2001; Ney and Graf, 2018; Sikma et al., 2020). It can be assumed that the sensible heat flux within the dense plant stand was largely determined by the entire canopy. In other words, during the day, mixing near the soil surface was impeded by stable temperature stratification, while in the evening cooling expanded upwards from the soil surface (Fig. 7f). In general, the processes described above were more pronounced during IOP 2 with its greater canopy height than with the lower canopy during IOP 1. The vertical wind profile showed consistently low wind speeds within the dense canopy (<0.5 m s⁻¹). Above the canopy layer, the wind speed increased in a log-like profile up to a maximum of 2 m s⁻¹.

3.2.2 Profiles of gross primary production

The detailed profile observations presented in the previous section enable us to calculate height-resolved estimates of gross primary production A. Using the 30 min averages of the vertical profiles for temperature, moisture, and CO₂ in the canopy, A is determined using the A-g_s model (Jacobs et al., 1997; Ronda et al., 2001). A (mg m⁻² s⁻¹) is calculated as follows:

$$\begin{array}{}\text{(2)}& A=\mathrm{LAD}\left(A_{\mathrm{m}}\right(h)+R_{\mathrm{d}}(h\left)\right)\left[\mathrm{1}-\exp \left({\displaystyle \frac{-\mathit{\alpha }\mathrm{PAR}\left(h\right)}{A_{\mathrm{m}}\left(h\right)+R_{\mathrm{d}}\left(h\right)}}\right)\right],\end{array}$$

where LAD (m${}_{\mathrm{leaf}}^{\mathrm{2}}$ m⁻³) is the leaf area density, A_m(h) is the CO₂ primary productivity (mg m${}_{\mathrm{leaf}}^{\mathrm{2}}$ s⁻¹) as a function of height h, R_d(h) (mg m${}_{\mathrm{leaf}}^{\mathrm{2}}$ s⁻¹) is the CO₂ dark respiration as a function of h, α (mg J⁻¹) is the light use efficiency, and PAR(h) (W m${}_{\mathrm{leaf}}^{-\mathrm{2}}$) is the amount of available photosynthetically active radiation within the canopy. Solar-zenith-angle-related variation in PAR intrusion and differences between atmospheric and skin values for temperature, moisture, and CO₂ are neglected. Figure 8a shows the winter wheat LAD applied in the calculation.

Figure 8 (a) Leaf area density (m${}_{\mathrm{leaf}}^{\mathrm{2}}$ m⁻³) on 7 May 2018 as a function of height (h) normalized to the maximum canopy h (h_c). The profile is typical for winter wheat, as defined by Olesen et al. (2004). (b) Time evolution of CO₂ gross primary production A (mg m⁻³ s⁻¹) on 7 May 2018 as function of h normalized to h_c. The profile is obtained using the profile measurements and using Eq. (2).

Figure 8b shows that the entire canopy contributes to the photosynthetic activity but with maximum A at $h/h_{\mathrm{c}}=\mathrm{0.7}$ (h_c: canopy height). This is primarily caused by the extinction of PAR within the canopy and reduced leaf density distribution close to the ground (Fig. 8a). Maximum diurnal productivity is found at around $h/h_{\mathrm{c}}=\mathrm{0.7}$, with the diurnal maximum at 12:00 UTC. Integration over the canopy shows minor discrepancies with respect to the bulk A-g_s model calculation, as the profile data allows for a more precise evaluation of photosynthetic activity. The profile measurements combined with Eq. (2) therefore allows for an improved modelling of the photosynthetic CO₂ uptake of vegetation depending on height and the understanding of mechanisms. More accurate estimates of CO₂ gross primary production still require improved knowledge of plant canopy micrometeorology (Drewry et al., 2014).

Figure 9 Source partitioning results for (a) H₂O and (b) CO₂ fluxes for IOP 1. Dashed grey lines show the measured latent heat flux (ET_ec) and net ecosystem exchange (NEE_ec) in half-hourly time steps. Values with subscript index p indicate estimates based on inversed profile concentration measurements (Sect. 3.4). Error bars for evaporation calculated from microlysimeters (E_lysi) and soil respiration measurements (R_s,ch) indicated to 1 SD. (ET_ec: evapotranspiration measured as latent heat flux L_vE by the eddy covariance system; E: evaporation; Tr_p: transpiration; NPP: above-ground net primary production; R_s: soil respiration).

3.2.3 Profile-based partitioning of H₂O and CO₂

Figure 9 shows the measured fluxes of latent heat, NEE, and soil respiration, as well as their partitioning based on the inversion of vertical high-resolution concentration profiles into the soil evaporation and plant transpiration and the R_s and NPP components. In this section, positive values indicate a flux from the surface and plants into the atmosphere and vice versa. During IOP 1, measured latent heat flux (L_vE, hereafter referred to as ET_ec) showed a typical daily pattern under clear sky conditions (Fig. 9a) with maximum ET_ec at noon (345 W m⁻²). Evaporation E of both methods displayed comparable values in the morning and evening but differed at midday. In the morning, the evaporation estimated using the profile measurements and method (E_p) and the lysimeter observations (E_lysi) both consistently suggested low E∕ ET fractions with E below 10 W m⁻². Towards noon, E_p increased to 25 and E_lysi to 60 W m⁻², and in the afternoon E_lysi reached a maximum of 101±41 W m⁻² (no E_p available). Estimated Tr_p increased to about 290 Wm⁻² at 11:00 UTC, this being the highest diurnal proportion of ET. Lower Tr_p levels around 12:00 UTC are probably due to a sub-optimal performance of the profile-based partitioning at this particular time. For example, none of the available inversion methods, including the algorithm by Santos et al. (2011) used here, include the effect of local thermal stability varying with height. Figure 7 demonstrates that thermal stability increased from the canopy top towards the ground around noon of IOP 1 (Fig. 7e), which may have contributed to the large increase of humidity towards the surface (Fig. 7b) due to the lack of mixing.

Figure 10 (a) Net primary production (NPP) versus photosynthetic active radiation (PAR). (b) Latent heat (L_vE) versus net radiation (Q^*). In both figures, the observation period encompasses clear and cloudy skies during a 2-week period starting on 7 May 2018 at 03:30 UTC (sunrise) and ending on the 20 May 2018 at 19:40 UTC (sunset). The solid line represents the median of the data. The lower and upper boundaries of the shaded area are the 25th and 75th percentiles of the data, respectively.

Variations in CO₂ fluxes NEE, NPP, and R_s during IOP 1 are shown in Fig. 9b. NEE_ec followed a typical diurnal cycle, with strong negative fluxes during the day and slightly positive values (carbon source) during transition times. The highest NEE was observed before noon (−25 µmol m⁻² s⁻¹). NPP_p followed the graph of NEE_ec, with higher values (−26 µmol m⁻² s⁻¹) in the morning hours than during the afternoon under comparable PAR values. This behaviour coincides with the photosynthesis rate observed at leaf level in Fig. 6a and provides further evidence that carbon uptake by plants was limited due to stomatal occlusion caused by the increase in VPD (Fig. 6a) and/or T_air in the afternoon. Profile-based R_s,p ranged between 0.5 and 6 µmol m⁻² s⁻¹ with higher values around noon. Compared to measured R_s,ch, R_s,p lay within the standard deviations of R_s,ch, though R_s,p was significantly lower during the morning and evening hours.

3.3 Effects of clouds on surface turbulent fluxes

3.3.1 Cloud-induced diffuse fertilization effect on evapotranspiration

One of the main aims of CloudRoots was to obtain observational evidence of the effects of clouds on the CO₂ assimilation and ET. Figure 10 shows the net primary production (NPP) (Fig. 10a) and L_vE (Fig. 10b), both measured using the eddy covariance, observed under a wide range of clear and cloudy skies as a function of PAR and compared to Q^* at the top of the canopy (van Diepen and Moene, 2019). We analyse a 2-week period of observations, between 7 and 20 May 2018. The effect of the different direct and diffuse radiation due to cloud perturbations is distinguishable with an enhancement of NPP under clear conditions whereas L_vE is reduced. Clouds affect plant photosynthesis by increasing the fraction of diffuse solar radiation that arrives at the top of the canopy (Kanniah et al., 2012). With a larger contribution of diffuse solar radiation and within the canopy, the radiation spreads more equally over all leaves and thereby increases the light-use efficiency of a canopy (Farquhar and Roderick, 2003). At a constant level of radiation at the top of the canopy, the increased light-use efficiency results in enhanced canopy photosynthesis, which is known as the diffuse fertilization effect (Roderick et al., 2001). This phenomenon is especially noticeable for canopies with a high LAI (Knohl and Baldocchi, 2008; Dengel and Grace, 2010). In CloudRoots and due to the high values of LAI (values between 4.5 and 5.5 m⁻²), we expect situations in which diffuse fertilization occurs, but here the question is how it influences L_vE. Previous large-eddy simulation modelling studies by Pedruzo-Bagazgoitia et al. (2017) have shown that under conditions dominated by clouds with a small optical depth, i.e. thin clouds, L_vE is enhanced with respect to its clear-sky values at the same radiation level.

Figure 11 IOP 2 (15 June 2018) time series of (a) sensible heat flux (H) at 1 min intervals with a displaced-beam laser scintillometer (DBLS) and at 10 min intervals with an eddy covariance system (EC), combined with scaled time series of photosynthetically active radiation (PAR, scaled by 1500 µmol m⁻² s⁻¹) and wind speed (U, scaled by 6 ms⁻¹); (b) turbulent exchange coefficient K_T; and (c) structure parameter of temperature, $C_T^{\mathrm{2}}$, that together make up H in the DBLS method following Eq. (1).

We find that the observed L_vE is higher rather than lower during clear conditions (less diffuse light) than under more diffused cloudy conditions. At constant Q^*, the median of L_vE is always higher under clear skies than for cloudy skies. The diffuse fraction plays a minor role and the decrease on L_vE under cloudy conditions is mainly due to the reduction in the incoming shortwave radiation. Our observations indicate that L_vE is driven by the partitioning of direct and diffuse radiation but also other effects such as diurnal variations of temperature, and the link to VPD may partially compensate for the different distribution of direct and diffuse radiation caused by clouds. The higher VPD values during the day partly offset the more optimal PAR conditions and therefore cause a closing of the stomatal that leads to decreases in L_vE. For both clear and cloudy skies, the shaded area below the median represents conditions before 11:30 UTC and the shaded area above the median represents conditions after 11:30 UTC, i.e. implying a hysteresis loop (Zhang et al., 2014). This spread in L_vE at a constant level of Q^* is caused by a difference in VPD between morning (before 11:30 UTC) and afternoon (after 11:30 UTC). This is because on a clear day the VPD raised rapidly due to its non-linear dependence on temperature relative to a cloudy day. In a typical clear day at CloudRoots, the value of 200 W m⁻² for Q^* is crossed twice: once in the morning and once in the afternoon. When 200 W m⁻² is crossed in the morning, the VPD is around 1000 Pa and reaches a value of 2000 Pa in the afternoon. On the other hand, on a cloudy day with similar values of around 200 W m⁻² the VPD remains almost constant through the entire day and with a value of 1000 Pa at 11:30 UTC.

The influence of VPD on L_vE also has the effect that the diurnal cycles of Q^* and L_vE are out of phase due to its dependence on leaf temperature. Q^* is primarily a function of incoming shortwave radiation and VPD of air temperature at the leaf surface. As a result, Q^* and VPD peak at different times of the day. Q^* peaks at maximum incoming shortwave radiation (local noon is at 11:30 UTC), and near-surface VPD times when air temperature peaks, which is around the time at which $Q^{*}=\mathrm{0}$ (17:00 UTC). The diurnal cycle of the sun implies there is a short period around 11:30 UTC when Q^* does not change. On the contrary, air temperature increases almost linearly around 11:30 UTC due to the approximately constant Q^*, as does VPD. Therefore, peak values for L_vE are found between the moments of maximum Q^* and of maximum VPD. For this dataset, the peak of L_vE is around 12:00 UTC for both clear and cloudy skies, although the peak for cloudy skies is less distinct due to the more fluctuating daily cycle of Q^*. Because Q^* and L_vE are out of phase, the highest values for L_vE do not occur in the bin with the highest net radiation, but rather in the bin of 400–500 W m⁻² (which roughly contains data from 11:00 UTC and after 12:00 UTC).

3.3.2 Cloud-induced radiation perturbations and response by turbulent fluxes

The short interval fluxes (1 min) of the double-beam laser scintillometer (DBLS) technique enable us to study the vegetation response to rapid radiation perturbations due to changes in cloud cover. The goal here is to illustrate this potential by discussing selected time series under changing cloud conditions during IOP 2. The morning of IOP 2 was characterized by rapidly changing cloud conditions due to the overpass of a shallow cumulus cloud deck. A breakdown of the 1 min DBLS sensible heat flux in terms of contributions from turbulent exchange (K_T) and the measure for temperature fluctuations ($C_T^{\mathrm{2}}$) is given in Fig. 11. This figure also depicts, on the same axes, scaled time series of wind speed and PAR that can be regarded as proxies that fuel mechanically induced turbulence (wind speed) and buoyancy turbulence (radiation in general) as well as photosynthesis (PAR).

Figure 12 IOP 2 (15 June 2018) time-series of (a) latent heat fluxes (LvE) at 1 min intervals with a displaced-beam laser scintillometer (DBLS) and at 10 min intervals with an eddy covariance system (EC) combined with scaled time series of photosynthetically active radiation (PAR, scaled by 1500 µmol m⁻² s⁻¹) and wind speed (U, scaled by 6 ms⁻¹). Panel (b) is the same as (a) but for the CO₂ flux ($F_{{\mathrm{CO}}_{\mathrm{2}}}$).

First of all, the 1 min DBLS fluxes of H closely follow the cloud cover induced radiation changes, but with a time-lag of 45–120 s (Fig. 11a). This is similar to those reported by van Kesteren et al. (2013b). H fluxes measured with EC techniques even when estimated over the relatively short interval of 10 min, which is not a standard output, are not capable of capturing such rapid dynamic behaviour of the flux regime (Fig. 11a). The dynamic behaviour in the DBLS H is mainly governed by fluctuations in T expressed by $C_T^{\mathrm{2}}$ (Fig. 11c) and to a lesser extent by changes in the exchange coefficient K_T (Fig. 11b). Note that is impossible to fully distinguish the three variables H, K_T, and $C_T^{\mathrm{2}}$ from each other as they are all interconnected, e.g. K_T is defined in terms of the Obukhov length L, which in turn depends on H and u_*. Nevertheless, our high-time-resolution observations demonstrate that changes in PAR induce very fast responses of the transported quantity T (Fig. 11c). Even in the absence of strong wind-induced variations in K_T, these T variations lead to approximately similar dynamic behaviour of H. On top of this, the additional but smaller wind-induced fluctuations in K_T are also reflected in H and lead to “noise” in the variability of H compared to the cloud-induced on–off behaviour of PAR.

Figure 13 (a) Aerial map of F₇₆₀ on 26 June 2018 with measurement locations used to combine with mobile (circles) and stationary (triangles) measurements. (b) Diurnal changes in F₇₆₀ on different days of the campaign as 5 min measurement averages depicted in the same colours as observation locations in (a).

Next we examine how soon the fluxes of H₂O and CO₂ respond to the cloud induced radiation changes. Figure 12 demonstrates that there is indeed a fast response, and the 1 min resolution fluxes of H₂O and CO₂ allow us to precisely determine a delay time of approximately 2 min for the increases in CO₂ uptake and transpiration of H₂O relative to the changes in PAR. The delay is once again undetectable with the standard 30 min eddy covariance results (Fig. 12). This behaviour is in line with what was concluded about the state of the vegetation observed at leaf level (Sect. 3.1). As the vegetation is not water-stressed and is at a stage of development at which it is still actively growing, it will react rapidly to changes in radiation, i.e. it is in a radiation-limited regime. Under the conditions of our study, stomata appear to have reacted only slowly or remained constantly open because leaves were unstressed or reacting only slowly to cloud-induced changes. Moreover, the timescale of a light-induced stomatal response (maximum values of 20 min; Van Kesteren, 2013b) is normally larger than the timescale of most fluctuations in radiation. Our suggested explanation is that the 1 to 2 min delay time observed between radiation and turbulent fluxes is due to processes associated with the inertia of the leaf in addition to turbulent transport between the leaf and laser path due to, e.g. the small but not negligible storage of heat, H₂O and CO₂ in the canopy layer. However, we need further evidence to disentangle the separation in delays between H₂O and CO₂ fluxes.

3.4 Sun-induced fluorescence (SIF) measurements: temporal variability

Studying spatial and seasonal variabilities in ET during plant growth was one of the key goals of CloudRoots. To this end, we analysed SIF observations measured on time and on space. The top-of-canopy measurements of SIF were carried out in two ways: (i) diurnal courses from a single representative location were recorded from a stationary FLOX system, and (ii) mobile measurements covering several locations within a field were recorded from a FLOX system that was housed in a backpack. To ensure reproducible measurements, the two fibre optics of the system were attached to a gimbal and were placed with a movable tripod 2 m above ground. Diurnal curves were acquired on 7 May, 4 June, and 14 June (only morning hours due to cloudy conditions in afternoon); mobile measurements (with change of measurement locations during the day) on 6 and 26 June. As SIF measurements should be performed under clear-sky conditions only, records affected by clouds were carefully removed. Aerial maps of SIF were acquired with the high-resolution imaging spectrometer HyPlant. Figure 13a shows the aerial map of F₇₆₀ acquired on 26 June, suggesting homogeneous canopy properties within the winter wheat study field, while great differences can be seen between different fields. The same image identifies the FloxBox measurement locations in the same colour code that reconstruct the diurnal temporal variability of F₇₆₀ during the entire CloudRoots campaign in Fig. 13b.

Diurnal changes in photosynthetic activity are clearly visible in F₇₆₀. Measurements made at different locations generally follow the same diurnal pattern, especially within the period 7 May to 14 June, further confirming the hypothesis that ET spatial heterogeneity within the winter wheat field was small. The seasonal changes are also traced by F₇₆₀: from 7 May until 14 June, the winter wheat canopy was photosynthetically active in a transition stage from booting (7 May) until grain filling (14 June), as is reflected by high SIF values. At the end of June, however, the canopy approached senescence and the reduction in photosynthesis was documented by greatly reduced fluorescence levels (see Fig. 13b, see pink values after 12:00 UTC). This photosynthesis reduction is also corroborated by the normalised difference vegetation index (NDVI), which was calculated as the normalized difference between far-red to red reflectance (see Supplement for details). The green dense canopy has a NDVI value close to 1, and the decrease in NDVI is caused by the yellowish colour of the winter wheat canopy (see Fig. S2 at the Supplement).

3.5 Connecting SIF and evapotranspiration flux at the landscape scale

It is difficult to directly quantify spatial variations in the ET flux with the currently available in-situ equipment due to the necessity of installing a large number of measurement stations. Recently, some promising concepts have been published that exploit the relationship between SIF and plant water relations (Damm et al., 2018; Jonard et al., 2020). Following these concepts, we studied the connections between ET to regional measurements of SIF in two steps, which were recorded on this scale by the airborne sensor HyPlant (see Fig. 13a). First, to obtain an estimation of the spatial variability ET at CloudRoots, we used the 15 km × 15 km map acquired by the HyPlant sensor on 26 June 2018 and a land use classification of the region (Lussem, 2018). ET cannot directly be measured, thus, it was predicted using different coefficients (K_c) that depend on the land use categories around CloudRoots. We define K_c as the ratio of ET over a particular crop relative to the ET of potential grass used as reference (Allen et al., 1998; Bogena et al., 2010). For this analysis, the regional land use map that consisted of 32 different land use classes was translated to a reduced classification scheme of 9 land use classes, which covered most of the vegetation types in the study region (Table 4). Roads were excluded from the analyses, as we assumed that their effect is negligible on the 15 m × 15 m grid.

Table 4 Estimated K_c coefficients for different land use classes that are dominant in the study area. The land use classes were calculated using a more detailed land use classification that consisted of 32 classes. For this study several classes that have similar transpiration rates were combined.

For the estimation of K_c coefficients to calculate ET, we used the plant developmental stage at the CloudRoots site at the end of June. For the main regional crops, namely sugar beet, winter wheat, winter barley, and potatoes, local measurements of evapotranspiration by EC towers were used. These data have been collected over several years and weekly averaged. This enabled us to compute K_c from measured and potential ET averaged over the last 2 weeks of June. In the particular cases of winter wheat and winter barley especially, the K_c coefficient changes rapidly at this time of the year, in extreme cases from 1.0 to 0.3 within 2 weeks, due to the onset of senescence. Therefore, the coefficients for these two crops shall be used with care. In absence of eddy covariance data, we calculates the characteristic values of K_c for each crop type and the developmental stage were taken from Allen et al. (1998). All estimated K_c coefficients for different crops can be found in Table 4. To estimate the ET over a specific area occupied by particular crop on a given day and time, the land use map was transferred to the map of K_c coefficients according to Table 4 and then multiplied by the potential ET, using the ET grass as a reference value (ETgrass), specific to that moment in time. Figure 14 shows the spatial variability of predicted ET for the IOP 3 inferred from the K_c coefficients and the value of potential grass reference averaged between 09:00 and 14:00 UTC. The area is a 1 km × 1 km square, characterized by a mean of 5.76 mmol m⁻² s⁻¹ and a standard deviation of 1.86 mmol m⁻² s⁻¹. Figure 14 shows that this method can provide plausible information on the variability of ET at the sub-kilometre scale, and it points out to the need to introduce this sub-grid ET variability information in modelling studies. In the second step of the procedure, we compared this estimated ET to the SIF measurements (F₇₆₀). Figure 15 shows the correlation between estimated ET and solar-induced fluorescence F₇₆₀ for 26 June (Julian day 177) for the different land covers. The correlation between mean F₇₆₀ values and predicted ET values is R²=0.61 with larger ET and F₇₆₀ vales for crops and grass compared to the forest conditions. It is calculated from the comparison pixel by pixel of the SIF (Fig. 13a) and ET (Fig. 14). As the HyPlant overflight was carried out at noon in order to acquire the maximal SIF values and minimize the influence of changing sun angle, we also used the maximal value of ET grass, measured at midday on 26 June. The large range of values of ET, F₇₆₀, and F₆₈₇ from the different land use categories corroborate the large variability of ET around the CloudRoots field.

Figure 14 Spatial variability of evapotranspiration inferred from combining K_c coefficients with the value of potential grass reference ETgrass. The x and y axes represent the geographical coordinates of the CloudRoots site in metres (50^∘52^′09^′′ N, 6^∘27^′01^′′ E).

3.6 Boundary layer integrated dynamics over heterogeneous landscapes

To integrate and improve the interpretation of our observations, we used CLASS to model the cloudless day 7 May 2018 (IOP 1). Our specific aims, related to the scales and processes under study, are (i) at leaf level to make use of the new constants in the mechanistic A-g_s model obtained from the observations (Fig. 5 and Table 3), (ii) at a landscape scale to represent the sensible heat flux in a heterogeneous landscape, and (iii) to estimate the potential impact of advection (heat) on the diurnal evolution of surface and boundary layer variables. Table A2 summarizes all initial and boundary conditions, constrained by the observations, which are employed in the modelling of the surface and atmospheric variables. Figure 16 compares the model results with the surface and upper-air observations. Focusing first on Fig. 16a, we found that the modelled H largely overestimates the observations taken at the CloudRoots. However, comparing our modelled H with the estimate of the regional flux shown in Fig. 1b, we found a satisfactory agreement in terms of magnitude and diurnal variability between this regional observed flux and CLASS model calculation. Note that here, and compared to Table 4, we oversimplified the land surface into two categories: “bare soil” and “vegetated”. To complete this evaluation, we show in Fig. S1 the impact of the optimized A-g_s constants presented in Table 3 (CloudRoots) versus the default ones. Both the evolution of surface fluxes and boundary layer height are in better agreement with the observations. Similar impacts on how leaf processes (rice) can influence the meteorology were reported by Ikawa et al. (2018). There the boundary layer temperature changed by up to 0.5 K depending on the constants used in the leaf photosynthesis model.

Figure 15 Relation between evapotranspiration (ET) and fluorescence F₇₆₀, including the standard deviation for the nine land use categories defined in Table 4. The data were collected on 26 June 2018.

Figure 16 Comparison of the model and observed results of 7 May 2018: (a) surface fluxes and (b) boundary layer depth. The regional H, an aggregate that combines the vegetated and bare soil surfaces around the CloudRoots site as shown in Fig. 1b, is also included. For the boundary layer depth estimations, we used three different observational techniques. The lidar and microwave (MW) techniques were located at the JOYCE site facility. Solid and dashed lines represent the model results of surface fluxes and boundary layer height with and without imposing the advection of heat, respectively (Table A2), for completing the information on initial and boundary conditions.

Our explanation of the improved comparison between the observations and the CLASS results using the aggregated sensible heat flux is the following: in a heterogeneous landscape such as the location of CloudRoots (Fig. 1a), each surface type contributes its own latent and sensible heat fluxes. It is the landscape aggregate of heat fluxes (named regionally and shown with triangles in Fig. 16a and introduced in Fig. 1b) and more specifically the sensible heat flux that governs the boundary layer evolution in terms of height, potential temperature, specific humidity, and atmospheric constituents. Only by using this higher H do we obtain satisfactory agreement with the observed boundary layer height evolution, which reaches its maximum values at around 1500 m in the afternoon (Fig. 8b). This further emphasizes that the H measured with the EC instrument during CloudRoots is only representative of the specific measurement site (leaf and canopy scales). The landscape average is an aggregate of values of H made up of the mosaic of surfaces as shown in Fig. 1. As a consequence, it is this composite H rather than a local value of H that is the main driver of the boundary layer development (boundary layer scales). With regard to ET, the model results are in good agreement with the local CloudRoots observations. This indicates the secondary and more local role played by ET in the dynamics of boundary layer development. For studies focusing on the regional values of ET, it will be necessary to calculate landscape-scale aggregate following the same procedures as H, while for studies at the leaf and canopy scales the local observations of ET are representative. Focusing now on Fig. 16b, we found a satisfactory agreement between the modelled boundary layer height and the three independent observations made with three different instruments. In Fig. 16b, it is interesting to note that the ABL height (75 km) inferred by the radio sounding measurement collected more than 100 km from of the Cloud Roots site has values similar to those collected by the lidar located within a radius of 5 km from the CloudRoots site. We attribute these similar values to a boundary layer that is characterized by being spatially homogeneous and with a similar temporal evolution on the larger regional scale.

CLASS, besides solving the diurnal variability of the boundary layer dynamics and the state variables, offers the possibility of adding a large-scale contribution that represents the advection of heat and/or moisture (see Vilà-Guerau de Arellano et al., 2015). We have performed a sensitivity analysis to determine the role played by heat advection for the surface fluxes and the boundary layer development. In the specific case that is modelled on 7 May, we relate this advection of heat or moisture to the diurnal evolution of H contrast between the measurement site and its adjacent fields, i.e. horizontal transport of heat, moisture, or momentum is driven by secondary circulations induced by the different thermal characteristics of the fields around the CloudRoots site (Fig. 1a). More specifically, we prescribe an advective heat contribution to represent the horizontal transport of heat due to the thermal variability of the surface conditions. This term follows an exponential function (Table A2 in Appendix A) with maximum positive values of advection equal to 0.9 K h⁻¹ at midday. This advective term is imposed only on the mixed layer and not on the free troposphere. Figure 16 shows how this advection of warm air to the CloudRoots site influences the boundary layer height. Starting with H, warm advection leads to higher mixed-layer temperatures that reduce the gradient between the temperature at the surface and the atmosphere, and thus reduce H. We find an opposite effect on ET. The increase in temperature by advection of warm air leads to an increased atmospheric demand, and therefore enhances ET. With regard to the boundary layer height, we might suppose that a drop in H would lead to a decrease of the boundary layer growth. However, the modelled boundary layer height displays the opposite behaviour. This is because the lower H is partly offset by a decrease in the thermal inversion at the interface between the boundary layer and the free troposphere. Lower values of the difference in θ_v between the free troposphere and the mixed layer enable boundary layer air parcels to be more easily transported into the free troposphere, resulting in faster growth of the boundary layer. This is because of the virtual potential temperature between the environment and the parcel is effectively reduced. The CLASS model results show that this process is more important than the decrease in H at the surface, and it allows the boundary layer to grow deeper than in the numerical experiment in which the warm advection is omitted. These numerical sensitivity experiment analyses enable us to quantify how non-local processes, in particular the effects of the regional average H and of warm advection, influence the observations at the measurement site.

Figure 17 Summary of midday evapotranspiration collected using different instrumental techniques during (a) IOP 1, (b) IOP 2, and (c) IOP 3. ET fluxes (left y axis) and g_sw (right y axis) reflect the period from 09:00 to 14:00 UTC. Box plots denote the variability in 30 min measurement intervals, except for the 1 min scintillometer measurements. The central mark of each box indicates the median, and the bottom and top edges of the box indicate the 25th and 75th percentiles, respectively. The acronyms are eddy covariance (EC), scintillometer (Scint) with 30 min and 1 min averages, ET inferred from the lift profiles (Lift), sap flow, ET calculated with the CLASS model, and ET inferred from the sun-induced fluorescence (SIF).

4 Discussion

CloudRoots offers an integrated methodology that combines field experiments across spatial scales (from leaf to landscape) closely linked to the modelling of the diurnal variability of the soil–plant–atmosphere continuum. To frame the discussion and link all our observations at the various scales and modelling efforts, we present in Fig. 17 all the different estimates of ET obtained during the three IOPs, averaged between 09:00 and 14:00 UTC in order to avoid the morning and afternoon transitions. Plotted alongside the ET estimates, we showed the leaf-level measurement of g_sw to indicate the control of vegetation on canopy-level ET. The four instrumental techniques are sap flow, eddy covariance (EC), scintillometer (averaged over 30 min and 1 min), ET inferred by the profile lift measurements, and ET infrared from the SIF observations. The ET modelled by CLASS is also included for IOP 1.

In comparing ET from the three IOPs, we find significant differences in magnitude from different techniques. In general, the highest values of ET are observed during IOP 1. The three IOPs were characterized by differences in the stages of growth, from very active vegetation to senescent vegetation, and influenced by a range of weather conditions: IOP 1 cloudless, IOP 2 scattered and thick clouds, and IOP 3 shallow cumuli. It is surprising that the decay in the vegetation activity as quantified by the measurements of leaf conductivity (Fig. 3d, e, f) is less evident in differentiating IOP 3 (senescent stage) from the more active vegetation at IOP 1 and 2. Furthermore, we observed, moving from IOP 1 to IOP 3, a much stronger decline in g_sw, suggesting that stomatal closure compensated for increased atmospheric moisture demand.

Several conclusions can be drawn from this intercomparison of ET observations using different techniques. Firstly, we might expect that the EC and scintillometer measurements, both with larger footprint and the inclusion of the soil evaporation contribution, show a net total ET that is similar to or higher than that obtained by the sap flow measurements. Secondly, we observed a far more pronounced response in declining g_sw compared to all ET measurements. These results point to the need to more accurately measure the leaf energy balance to take the penetration of radiation in the canopy under clear and cloudy conditions into account. This would also require a revision of scaling procedure from the leaf to the canopy level. Secondly, it is known that the EC flux measurements normally underestimate the sensible and latent heat fluxes because the EC flux measurements filter out the low frequencies (Foken et al., 2008; Gao et al., 2017). This underestimation is difficult to determine, but as a first guess and related to Fig. 17, the underestimation might range between 10 % and 15 %.

Although the contribution of soil evaporation is small compared to plant transpiration due to the high vegetation cover, we need to stress that EC and scintillometer observations are similar to or smaller than the ET observed or inferred from the other techniques (Fig. 17). This highlights the difficulty of estimating ET due to the need to include and quantify the contributions of the four fundamental processes: soil evaporation, upscaled leaf transpiration, evaporation related to the sap flow and the two non-local processes, entrainment of dry air, and horizontal advection of heat and moisture. Here, the modelling of ET, taking into account for and integrating all these processes, enables us to discriminate among these processes and calculate the budget of ET as a function of these local and non-local contributions. In that respect, the CLASS model is a tool capable of efficiently combining observations and model results that integrate surface and boundary layer dynamics. The averaged modelled ET is at the higher range of the ET observed estimations during IOP 1.

With respect to the differences between the 1 min and 30 min series measured by the scintillometer, their median is very similar in the three IOPs. However, differences become larger at smaller timescales due to the non-steadiness of ET under the presence of clouds. Here, the 1 min flux calculated from the scintillometer can capture the rapid and large fluctuations by clouds (Fig. 12) and the maximum values in particular. In order to obtain more definitive conclusions on how ET varies under cloud conditions, we need to analyse other situations characterized by different diurnal cloud cycles in more detail and systematically relate ET to key cloud characteristics such as the cloud optimal depth to determine how cloud thickness influences ET and the timescale of the cloud passage.

Regarding the quantification of the different processes contributing to ET, Fig. 9 illustrates the need to continue to test analytical techniques to identify the individual contributions of soil and plants to determine the diurnal ET budget. A possibly useful tracer would be the stable isotopic composition of water vapour and carbon dioxide (Lee et al., 2009; Griffis 2013) and, combined with isotope signals, for modelling the surface and boundary layer dynamics with the carbon and water exchanges. To further discriminate between soil and plant sources and sinks under unsteady conditions due to radiation and dynamic perturbations by cloud shading, these high-frequency stable isotope measurements should go beyond the typical average time of eddy covariance (30 min). As van Kesteren et al. (2013) showed and is further corroborated in this work, the scintillometer technique combined with high-frequency observations of H₂O and CO₂ enable us to quantify the response time of ET and CO₂ assimilation to these intermittent radiation fluctuations or cloud flecks (Kaiser et al., 2018).

Finally, the integration of all processes in the CLASS model shows the challenges in interpreting the measurements taken at the sub-kilometre scales and adequately representing the surface turbulent fluxes. Although the measurements indicate that the day selected for the modelling displayed a very homogeneous boundary layer depth over an area with a radius of 100 km², the sensible heat flux measured at the CloudRoots facility was not representative of it. Therefore, we recommend extending the number of stations by means of a multi-tower approach that would also include detailed observations of the soil and plant conditions. In addition to obtaining a more representative field sensible heat flux that is better related to the development of the boundary layer, a denser network of spatial observation stations is also necessary to estimate the role of hectometre-scale heterogeneity-induced circulations and their relationships with the local advection of heat and moisture more accurately (Mauder et al., 2010).

5 Conclusions

Our main findings, organized from the smaller to the larger scales observed and modelled, are summarized as follows:

At leaf scale, we find that stomatal conductance and gross primary production decrease in line with the increasing senescence of the plant. The tiller-level measurements of the sap flow are virtually constant throughout the growing period. Underlying causes need to be further investigated under controlled conditions. The successful modelling of the leaf stomatal conductance and the photosynthesis assimilations required the relevant constants used in the mechanistic model (A-g_s) in the field to be measured. Modelled leaf-level photosynthesis compares better with the measurements during the mature growing period than during senescence. For future field experiments, we recommend including leaf-level measurements in meteorological campaigns to improve calculations related to the water–carbon leaf and canopy exchanges.

At canopy scale, the high-frequency vertical profiles – measured in and above the canopy – of wind speed, potential temperature, specific humidity, and carbon dioxide prove to be very valuable in obtaining profiles of gross primary production in the canopy and as a function of height. By inverting these observed profiles, we obtain an estimate of the contributions of soils and plants to the net evapotranspiration and CO₂ ecosystem exchange. The validation against individual measurements of these components gives better results for the net ecosystem exchange than those for the net evapotranspiration. We argue that for evapotranspiration the dependence on temperature and water vapour deficit plays a more important role than for CO₂ assimilation, the latter being mainly controlled by the partitioning between direct and diffuse radiation.

Under cloud conditions, we show that the perturbation by clouds of direct and diffuse radiation create large fluctuations in evapotranspiration and the CO₂ assimilation with opposite signs for evapotranspiration and CO₂ exchange. A cloudy boundary layer reduces evapotranspiration, whereas it enhances plant assimilation of CO₂. The 1 min turbulent fluxes acquired by the scintillometer demonstrate the relevance of flux measurements observed at higher frequencies for improving quantification of the impact of clouds on the photosynthetically active radiation. With these fast turbulent fluxes, we quantify delays of the turbulent fluxes with respect to the photosynthetically active radiation. These delays are on the order of minutes. Comparing these 1 min flux estimates with the standard 30 min average measured with the eddy covariance technique, we find a lower median and a large increase in the variability of the net evapotranspiration. This information can be useful in determining the impact of rapid fluctuations driven by the impact of clouds on evapotranspiration and its impact on the closure of the surface energy balance.

At landscape and boundary layer integrated scales, the modelled sensible heat flux correlates better with the area-weighted average flux than the local flux estimates. The area-weighted flux integrates in a simple manner a composite of bare soil and vegetated surfaces at regional scale (kilometres). This aggregate regional flux is representative of an area that is larger than the CloudRoots site (100 m × 100 m). Therefore, a model setup that represents the boundary layer evolution well only needed to be informed by the area-weighted average of two main surface types, bare soil, and vegetated areas. The variations of ET due to surface heterogeneity were also measured and inferred from airborne sun-induced fluorescence observations. Our findings corroborate the large heterogeneity of ET at the sub-kilometre scales with values ranging from forest (about 2.5 mmol m⁻² s⁻¹) to late crops such as potato or sugar beet (8–10 mmol m⁻² s⁻¹) .

The comparison of all the ET measurements at the various scales show that there are still large differences in observing ET among the different observing techniques, the modelling of ET and their relation to stomatal aperture during the entire growing season. These ET observations do not show a clear pattern related to the scale at which they were measured.

The modelling and scale integration of this comprehensive observational dataset enables us to study the carbon and water exchange at leaf, canopy, and landscape levels. It also allows us to quantify how horizontal advection of heat within the mixed layer influences the surface fluxes and the growth of the atmospheric boundary layer. We show, for instance, that the horizontal advection of heat leads to deeper boundary layer depths. This numerical experiment thus paves the way to more complete modelling studies, for instance using large-eddy simulation numerical experiments, on how the surface and the overlaying atmosphere interact on sub-diurnal and sub-kilometre scales.

Appendix A: Construction of light and radiation response curves under clear and cloudy conditions

For the construction of the light and radiation response curves in Fig. 10, the data were divided into bins of PAR and Q^*. For Fig. 10a, we divided the data points in bins of incoming total PAR. Each bin covers a range of 50 Wm⁻², starting at 0–50 Wm⁻² and presumably ending at around 350–400 Wm⁻² (maximal intensity for incoming direct PAR). For Fig. 10b and the variable Q^*, each bin covers a range of 100 W m⁻², starting at −200 to 100 W m⁻² for cloudy skies and at 0 to 100 W m⁻² for clear skies. In both figures, for each data point the diffuse fraction of PAR is determined by combining measurements of incoming total PAR and incoming diffuse PAR. Subsequently, a data point is labelled “clear” for diffuse fractions <0.3 and labelled “cloudy” for diffuse fractions >0.8. We choose these boundaries to balance a distinct difference between clear and cloudy skies with a large enough sample size for each bin. For clear skies, the first two bins are missing. This is due to the fact that under clear skies low levels of Q^* are the result of the sun being close to the horizon, and as a result solar radiation has to travel a long distance through the atmosphere before reaching the surface. In those cases, most of the solar radiation reaches the surface as diffuse radiation due to Rayleigh scattering and scattering by aerosols, and therefore does not meet the criteria to be labelled “clear”. For cloudy skies, bins are missing for high levels of Q^*. Clouds attenuate solar radiation through absorption and backscattering, and thereby reduce Q^* to a level lower than it would be for a clear sky.

Table A1 Management activities on the test site over the winter wheat cultivation cycle before, during, and after the observation period of the CloudRoots campaign.

Table A2 Initial and boundary conditions prescribed in CLASS to reproduce IOP 1 (7 May 2018).