Surface flux estimates derived from UAS-based mole fraction measurements by means of a nocturnal boundary layer budget approach

The carbon exchange between ecosystems and the atmosphere has a large influence on the Earth system and specifically on the climate. This exchange is therefore being studied intensively, often using the eddy covariance (EC) technique. EC measurements provide reliable results under turbulent atmospheric conditions, but under calm and stable conditions – as they often occur at night – these measurements are known to misrepresent exchange fluxes. Nocturnal boundary layer (NBL) budgets can provide independent flux estimates under stable conditions, but their application so far has been limited by rather high cost and practical difficulties. Unmanned aircraft systems (UASs) equipped with trace gas analysers have the potential to make this method more accessible. We present the methodology and results of a proof-of-concept study carried out during the ScaleX 2016 campaign. Successive vertical profiles of carbon dioxide dry-air mole fraction in the NBL were taken with a compact analyser carried by a UAS. We estimate an average carbon dioxide flux of 12 $\mathrm{\mu }\mathrm{mol}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{s}}^{-\mathrm{1}}$, which is plausible for nocturnal respiration in this region in summer. Transport modelling suggests that the NBL budgets represent an area on the order of 100 km².

1 Introduction

The terrestrial biosphere plays a major role in the carbon cycle. It has taken up approximately one-quarter of the carbon emitted by human activities since 1750 (Ciais et al., 2013), but the future development of this land sink under a changing climate is uncertain. Given its importance, the biosphere–atmosphere exchange is being studied intensively. On the ecosystem level, sources and sinks of carbon dioxide are commonly quantified using the eddy covariance (EC) technique (Baldocchi, 2003). During the day, when the air is turbulently mixed, the EC technique provides reliable direct measurements of net ecosystem exchange (NEE). However, EC measurements often misrepresent nighttime fluxes (Goulden et al., 1996; Gu et al., 2005). This is related to the stable stratification that often develops close to the surface at night. Stable conditions combined with low wind speeds violate assumptions underlying the EC technique (see Aubinet et al., 2012, for a comprehensive discussion). Despite large efforts, there is currently no generally accepted solution for how to obtain reliable measurements of nighttime fluxes using the EC technique (Gu et al., 2005; Aubinet et al., 2010; Hayek et al., 2018).

Daytime NEE consists of photosynthetic uptake and release of carbon through respiration. Nighttime NEE is governed by respiration only, as photosynthesis cannot take place without light. Photosynthetic uptake and total respiration fluxes are usually of the same order of magnitude, but with opposite sign. Therefore, even slight underestimation of nocturnal respiration can result in a considerable overestimation of an ecosystem's long-term carbon uptake. Furthermore, daytime fluxes are often partitioned into photosynthetic uptake and respiration using methods that rely on the nighttime measurements (Falge et al., 2001; Reichstein et al., 2005; Lasslop et al., 2010; Wohlfahrt and Galvagno, 2017). Errors in the nocturnal fluxes might compromise this partitioning.

The nighttime problem of EC measurements calls for error quantification and potentially correction. Ideally, this would be achieved by comparison to a method that is sensitive to fluxes in the same area but is not based on the same assumptions as the EC technique. Unfortunately, no such method is available. Constraining the error of EC measurements must therefore rely on methods that determine fluxes on smaller and larger scales.

Enclosure-based methods and biometric approaches, including plant growth assessment and stock inventories, are often employed to obtain independent estimates for NEE (Goulden et al., 1996; Wilson and Baldocchi, 2001; Campioli et al., 2016). These methods quantify the exchange of carbon on a much smaller spatial scale than EC measurements. The chambers typically used for determining soil respiration cover an area of less than 1 m², while the EC technique is sensitive to fluxes from an area of 10⁴–10⁶ m², depending on the site and on meteorological conditions (Chen et al., 2009). Given these different scales, inhomogeneities in the ecosystem under study, such as spatial variability of soil properties (e.g. texture, carbon content, nitrogen content), soil environmental conditions (e.g. soil temperature and moisture) or plant community composition, can lead to biases in the comparison.

In order to counteract these biases, larger-scale flux estimates should be obtained in addition to enclosure-based or biometric measurements when constraining the error of EC measurements. Nocturnal boundary layer (NBL) budgets, first described by Choularton et al. (1995) and Denmead et al. (1996), provide such estimates. The NBL budget method makes use of the stable stratification at night, which can act as a flux-integrating enclosure. During clear nights, the emission of thermal radiation cools down the Earth's surface much faster than the air, owing to the surface's higher emissivity. An inversion layer forms, inhibiting exchange of air between the stable NBL and the neutral residual layer above (Stull, 1988). Any tracer emitted from the surface into the atmosphere is therefore accumulated within the NBL. By measuring the rate of accumulation, the tracer flux can be estimated.

Different setups have been used for NBL budgeting. Acevedo et al. (2004) measured CO₂ dry-air mole fractions at a 12 m tower, sampling only the lowest parts of the NBL. Although they did not sample the whole layer, they were able to determine a budget by identifying an effective accumulation height from either heat flux or balloon-borne humidity and temperature measurements and assuming a uniform accumulation rate of CO₂ up to this height. Winderlich et al. (2014) used CO₂ and CH₄ dry-air mole fraction measurements at six heights on a 301 m tall tower, yielding profiles that encompass the whole NBL during most nights. Often the NBL budget method is applied without a tower. A tethered balloon can lift a 100–300 m long hose through which a ground-based gas analyser samples air from different heights (Choularton et al., 1995; Denmead et al., 1996; Culf et al., 1999). Alternatively, a light analyser can be carried by the tethered balloon directly (Pattey et al., 2002).

Despite providing unique information, the NBL budget method has been applied only infrequently in recent years. This might be related to the cost and operational limits of towers and tethered balloons. Unmanned aircraft systems (UASs) could make the NBL budget method more accessible. UASs with payload capacities on the order of 1 kg are now available for a few thousand euros. When equipped with lightweight trace gas analysers and meteorological sensors (Kunz et al., 2018), they have the potential to probe the NBL with large flexibility at low cost. Multicopters are a particularly attractive type of UAS for this kind of study, because their vertical take-off, vertical landing and hovering capability makes them easy to operate in a range of environments.

However, the air movement caused by a UAS can disturb the NBL and thereby compromise the measurements. A reliable NBL budget can be determined only if this issue is addressed. A second challenge is not specific to UASs but common to all NBL budgets: the area contributing to the budget depends heavily on the meteorological conditions and can extend far from the point of measurement. For a given time and site this footprint cannot be influenced by experimental design. Nevertheless, knowledge of the footprint is beneficial for the data analysis and interpretation of the results. In earlier NBL studies, this topic received only basic treatment (Culf et al., 1999) or was ignored altogether.

To assess the suitability of UASs as measurement platforms for the creation of NBL budgets, we carried out a proof-of-concept study. We deployed a carbon dioxide analyser on a multicopter and repeatedly sampled vertical profiles of the NBL during two nights in July 2016 as part of the ScaleX 2016 campaign in Fendt, Germany. Section 2 is a description of the site and the available ground-based instrumentation, the airborne measurement system and the unmanned aircraft. In Sect. 3 we explain how we dealt with the disturbance caused by the UAS, which post-processing steps we carried out and how we determined the NBL budget. Furthermore, we delineate how a Lagrangian transport model can be applied to identify the areas that contributed to the budget, i.e. how to determine the footprint of our flux estimates. In Sect. 4 we present and discuss the profiles taken by the UAS, the fluxes obtained from the NBL budget and a summary of the footprint analysis. We compare our observations to references and assess the robustness of our flux calculation. In Sect. 5 we summarise the merits and experimental challenges of our approach.

2 Site and instrumentation

2.1 Fendt site

The Fendt site is located in southern Germany in the Alpine Foreland (Fig. 1) at 47.833^∘ N, 11.060^∘ E (WGS84), 600 m above mean sea level.

Figure 1 Location, topography and aerial image of the Fendt site and its surroundings. MOHP is the Meteorological Observatory Hohenpeißenberg, and HPB is the ICOS station Hohenpeißenberg. Digital elevation model and aerial imagery by Bayerische Vermessungsverwaltung, https://www.geodaten.bayern.de (last access: 25 March 2020).

The site lies in a flat valley bordered by a gentle slope to the east and a steep slope leading to a 100 m higher plateau to the west. The valley floor is dominated by pasture and some crops, predominantly maize, which in Germany is typically sown in April or May and harvested between September and November. The slopes to the east and west are covered with coniferous and mixed forest, respectively. Fendt belongs to the district Weilheim–Schongau, which has a population density of 139 km⁻² (Statistisches Bundesamt (Destatis), 2018, 35th percentile of all districts in Germany).

While soil identification at the Fendt site resulted in Stagnosols at three locations, soil organic carbon (SOC) content was determined additionally at 20 locations within a regular grid covering an area of 300 m by 300 m. SOC content at 5 cm depth varied between 4 % and 11 %, while at 50 cm depth, values of up to 23 % were obtained. The highest SOC contents were observed at the eastern side of the regular grid where a peat area is located. According to BGR (2013), organically rich soils (Cambisols and Histosols) prevail within a 20 km radius around the Fendt site (Fig. 2a). The dominant land cover in this region includes crops, pasture and forest (Fig. 2b).

Figure 2 (a) Soil types in the region around the Fendt site, based on BGR (2013), denoted in WRB classification (IUSS Working Group WRB, 2015). “Fl./Gl.” stands for “Fluvisols/Gleysols”. (b) Simplified land cover map (CORINE 2012 v18.5, European Environment Agency, EEA, 2016) of the same region. In both panels the location of the Fendt site is marked with a black diamond.

About 5 km south-west of the Fendt site lies an isolated, 988 m high mountain, the “Hoher Peißenberg”. Close to its summit the German Weather Service (Deutscher Wetterdienst, DWD) operates the Meteorological Observatory Hohenpeißenberg (MOHP) and the ICOS (Integrated Carbon Observation System) station Hohenpeißenberg (HPB, Fig. 1).

2.2 Ground-based instrumentation

Fendt is part of the TERENO (Terrestrial Environmental Observatories) network and is extensively instrumented for the purpose of long-term monitoring of land–atmosphere exchange (Mauder et al., 2013; Kiese et al., 2018). Complementary observations were made during the ScaleX 2016 campaign (Wolf et al., 2017). In the following, we list only those instruments that produced the data presented in this publication.

During the ScaleX 2016 campaign, CO₂ dry-air mole fraction at heights of 1, 3 and 9 m above ground level was measured with a cavity ring-down spectrometer by successive sampling of air through three inlets installed at a 9 m high mast. Each inlet was sampled once every 7.5 min, with occasional interruptions due to calibrations and other measurements. An EC system installed at 3.5 m height (Zeeman et al., 2017) quantified the turbulent exchange of CO₂. Air temperature as well as upward and downward radiation were measured at 2 m height. Two sets of automated chambers were operated to determine the total NEE or respiration flux of grass and soil. One set comprised four LI-8100 long-term chambers (LI-COR, Lincoln, NE, USA), two with a clear enclosure for measuring NEE and two with an opaque enclosure for measuring respiration (Zhao et al., 2018). All four chambers covered an area of 317.8 cm² and will be referred to as “small chambers” from here on. The other set consisted of five custom-built opaque chambers covering an area of 2500 cm², referred to as “big chambers” hereafter. All the instruments mentioned so far were located close to each other, and the horizontal distance between any of the instruments and our UAS was always smaller than 200 m during the flights.

Besides the on-site instruments, we use two more data sources for our analysis. One is the observation of cloudiness at the Meteorological Observatory Hohenpeißenberg, recorded every hour either by a person or an automated instrument. We consider these 5 km distant measurements representative for Fendt, with a potential time lag on the order of 1 h in the case of synoptic events. The second non-local data source is the greenhouse gas monitoring system at the ICOS station Hohenpeißenberg, situated at 934 m above mean sea level. We use its measurements of the CO₂ dry-air mole fraction at 131 m height above ground level, i.e. at 460 m above the Fendt site.

2.3 Airborne payload

For the study presented here, temperature, pressure, relative humidity and CO₂ dry-air mole fraction of ambient air were measured using COCAP, the COmpact Carbon dioxide analyser for Airborne Platforms, developed at the Max Planck Institute for Biogeochemistry in Jena (see Kunz et al., 2018, for a detailed description). COCAP was mounted below the multicopter. Air samples for the measurement of carbon dioxide dry-air mole fraction were drawn from an inlet placed 30 cm below and 20 cm to the side of the rotors (Fig. 3).

Figure 3 COCAP was carried by a multicopter during the ScaleX 2016 campaign. The positions of the sample inlet for the CO₂ measurement and of the temperature and humidity sensor board are indicated.

The temperature and humidity sensor board, requiring strong ventilation for the fastest response, was placed directly below one of the rotors. The sensor for ambient pressure was located inside COCAP's housing, which was not hermetically sealed and therefore in equilibrium with ambient pressure.

The measurement principles employed by the different sensors as well as their measurement uncertainties are listed in Table 1. The uncertainty of the calibration is included in the measurement uncertainties reported.

Table 1 Measurement principles, uncertainties and calibration range of the airborne sensors.

2.4 Unmanned aircraft

During ScaleX 2016 COCAP was deployed on an S1000 multicopter (SZ DJI Technology, China) controlled by a Pixhawk autopilot (3D Robotics, Berkeley, CA, USA) running the Ardupilot APM:Copter v.3.3.3 firmware. Take-off mass of the whole system was 8 kg. The multicopter was powered by three lithium polymer batteries with a voltage of 22.2 V and a capacity of 5000 mAh each, achieving a maximum flight time of 12 min. Our special flight permit included nighttime flights, but because the take-off mass of our UAS exceeded 5 kg, all flights were limited to a maximum height of 150 m.

3 Methods

3.1 Disturbance by the UAS

A multicopter as a rotary-wing aircraft counterbalances gravity by accelerating air downwards through the movement of its rotors. The resulting displacement of air can interfere with in situ measurements, because air might be sampled at a location where it would normally not reside. In addition, volumes of air originating from different locations can be mixed together. The greater the displacement and mixing caused by the UAS, the greater the potential impact on the measurement of a gradient, for example. Air movement below and above the rotors is not symmetric: below a rotor, air is pushed downwards as a directed stream with high speed. In contrast, the air flow towards the rotor comes from different directions and has a lower speed. The reader can easily confirm this with a fan or a hair dryer: while the outflow of air can be felt metres away, the inflow is hard to sense even near the rotor.

In view of the asymmetric flow pattern, we expect that during ascent of the UAS air parcels are measured with negligible displacement from their undisturbed location. During descent, however, the sensors are moved into a volume that potentially has been flushed with air originating from several metres above. During hovering at a fixed location or during purely horizontal movement, the sensors might reside in a partially closed flow loop that extends below and to the side of the multicopter, effectively measuring a mixture of air from different locations.

For the study presented here, flying near the ground can have a particularly strong influence on the measurements for three reasons. Firstly, downward motion of the air stops at the ground and displaced air must move laterally or upwards, making a fast flow path back to the UAS more likely. Secondly, in our nighttime experiments the air near the ground is stably stratified. Therefore, air pushed downwards by the rotors experiences a restoring upward force, increasing the chance that closed flow loops form. Thirdly, the strongest gradients in temperature and CO₂ dry-air mole fraction are present close to the ground, hence even a small displacement of air can have a large effect on the measured values.

In the case of considerable horizontal air speed, due to either wind or horizontal flight, the rotor-induced airflow should have a smaller effect on measurements because the sampling system is moving away from air that has been displaced. We investigated this effect by flying horizontally at different speeds over a homogeneous meadow (see Sect. 4.4).

Based on the considerations above and the data presented in Sect. 4.3 and 4.4, we determine the NBL budget only from those measurements that were taken during ascent of the multicopter. The sensitivity of the NBL-derived fluxes to inclusion of hover and descent data is discussed in Sect. 4.5. Furthermore we discard COCAP's $x_{{\mathrm{CO}}_{\mathrm{2}}}$ data collected below 9 m height for the calculation of the NBL budget. Instead, the lowest part of the $x_{{\mathrm{CO}}_{\mathrm{2}}}$ profile is defined by the stationary measurements at the 9 m mast at 1, 3 and 9 m height. Pressure and temperature at these levels are interpolated from COCAP's measurements. During flight, the horizontal distance between COCAP and the 9 m mast was less than 150 m at any time. Hence, we do not expect pronounced horizontal gradients in $x_{{\mathrm{CO}}_{\mathrm{2}}}$ between the measurement locations. In Sect. 4.6 we discuss how the NBL-derived fluxes are affected if the data from the 9 m mast are not used.

3.2 Correction for response time of sensors

On a moving platform the finite response time of sensors can be a source of measurement error, as the response time distorts the attribution of data points to time and location. COCAP's pressure and temperature sensors are fast enough for this effect to be neglected, but both the humidity and the CO₂ sensor require correction.

The response of a capacitive humidity sensor can be expressed following Miloshevich et al. (2004):

$$\begin{array}{}\text{(1)}& {\displaystyle \frac{\mathrm{d}{U}_{\mathrm{m}}}{\mathrm{d}t}}=k(U_{\mathrm{a}}-U_{\mathrm{m}}).\end{array}$$

U_a and U_m are ambient and measured relative humidity, respectively. The coefficient k is inversely related to the sensor's response time and might be temperature-dependent. Solving Eq. (1) for U_a provides a simple way to compute true humidity from measurements. We use a fourth-order Savitzky–Golay filter (Savitzky and Golay, 1964) with a length of 15 samples to compute $\frac{\mathrm{d}{U}_{\mathrm{m}}}{\mathrm{d}t}$ while keeping high-frequency noise at an acceptable level. The coefficient k=14 s was determined by an optimisation that minimises the difference between the corrected humidity profiles for ascent and descent. We tested a linear and quadratic dependence of k on ambient temperature, but we found no improvement that would justify the additional degrees of freedom in the model.

The response of COCAP's CO₂ sensor is more complex. Its response to step changes in CO₂ dry-air mole fraction can be approximated as

$$\begin{array}{}\text{(2)}& x_{\mathrm{SC}}\left(t\right)=\left\{\begin{array}{ll}{x}_{\mathrm{0}}& \text{if }t<t_{\mathrm{d}}\\ a(x_{\mathrm{0}}-x_{\mathrm{\infty }})e^{(t_{\mathrm{d}}-t)/{\mathit{\tau }}_{\mathrm{1}}}+(\mathrm{1}-a)& \\ \cdot (x_{\mathrm{0}}-x_{\mathrm{\infty }})e^{(t_{\mathrm{d}}-t)/{\mathit{\tau }}_{\mathrm{2}}}+x_{\mathrm{\infty }}& \text{if }t\ge t_{\mathrm{d}}.\end{array}\right.\end{array}$$

Here x₀ and x_∞ denote the CO₂ dry-air mole fraction before and infinitely long after the step change, respectively, and t_d is the sensor's dead time.

We determined the coefficient a, the dead time t_d, and the time constants τ₁ and τ₂ from experimental data collected in the field. With COCAP running in flight configuration, i.e. with the inlet tube attached, we connected a tube with gas flowing from a cylinder. We observed a dead time of t_d=5 s between making the connection and the first change of COCAP's reading. The remaining parameters were found by least-squares regression of Eq. (2) to the data, yielding τ₁=27 s and τ₂=3.2 s.

Ignoring noise and calibration error, any CO₂ signal x_a is reported by COCAP as the convolution of x_a with the CO₂ sensor's instrument function f (see Kunz et al., 2018):

$$\begin{array}{}\text{(3)}& x_{\mathrm{m}}\left(t\right)=(x_{\mathrm{a}}\cdot f)\left(t\right)\end{array}$$

with

$$\begin{array}{}\text{(4)}& f=\left\{\begin{array}{ll}\mathrm{0}& \text{if }t<t_{\mathrm{d}}\\ \frac{\mathrm{1}}{x_{\mathrm{\infty }}-x_{\mathrm{0}}}\cdot \frac{\mathrm{d}{x}_{\mathrm{SC}}}{\mathrm{d}t}& \text{if }t\ge t_{\mathrm{d}}.\end{array}\right.\end{array}$$

As x_SC is known from experiments, f can be calculated. The ambient signal x_a can be recovered from the measured signal x_m by deconvolution (Fig. 4). We carried out the deconvolution in Fourier space where it is equal to a division. In the numerical implementation it is important to discretise f in a way that does not underestimate the slope of f between the time steps t_d and t_d+Δt, because doing so would lead to a strong enhancement of the noise during the deconvolution of x_m with f. The opposite error, i.e. overestimating the slope between t_d and t_d+Δt, is less critical and just results in a slight smoothing.

Figure 4 Response of COCAP to an abrupt change in CO₂ dry-air mole fraction $x_{{\mathrm{CO}}_{\mathrm{2}}}$ at the inlet at time t=0 s. The measured signal reveals a dead time of 5 s of the sampling system. Furthermore, the step change in $x_{{\mathrm{CO}}_{\mathrm{2}}}$ is smoothed out. Both effects are removed by deconvolution at the cost of higher noise. Smoothing the deconvoluted signal reduces the noise with only minor impact on the time response. Smoothing was carried out by convolution with a Gaussian function of 5 s full width at half maximum (FWHM).

3.3 Calculation of the NBL budget

For a parcel of air in the atmosphere the following continuity equation holds (Leuning, 2004, Eq. 6.2):

$$\begin{array}{}\text{(5)}& {\displaystyle }\mathit{\sigma }& {\displaystyle }={\displaystyle \frac{\partial c}{\partial t}}+\mathrm{\nabla }\cdot \left(c\mathit{u}\right),\text{(6)}& {\displaystyle }& {\displaystyle }={\displaystyle \frac{\partial c}{\partial t}}+c(\mathrm{\nabla }\cdot \mathit{u})+u{\displaystyle \frac{\partial c}{\partial x}}+v{\displaystyle \frac{\partial c}{\partial y}}+w{\displaystyle \frac{\partial c}{\partial z}}.\end{array}$$

Here σ is the strength of a volume source (or sink) of carbon dioxide ($\mathrm{\mu }\mathrm{mol}{\mathrm{m}}^{-\mathrm{3}}{\mathrm{s}}^{-\mathrm{1}}$), c is the concentration of carbon dioxide and t denotes time. The components u, v and w of the wind vector u point towards east (x direction¹), towards north (y direction) and upwards (z direction), respectively. Molecular diffusion is neglected. Due to continuity of the air flow, the term c(∇⋅u) equals zero. If we follow an air parcel as it is transported by horizontal winds, those terms that contain a horizontal wind component vanish as well and Eq. (6) is reduced to

$$\begin{array}{}\text{(7)}& \mathit{\sigma }={\displaystyle \frac{\partial c}{\partial t}}+w{\displaystyle \frac{\partial c}{\partial z}}.\end{array}$$

Now we integrate vertically over those air parcels that form a vertical column over our site at the time of measurement (at earlier or later points in time, the air parcels are not aligned in a vertical column, unless the wind vector is equal at all heights):

$$\begin{array}{}\text{(8)}& \underset{\mathrm{0}}{\overset{z_t}{\int }}\mathit{\sigma }\mathrm{d}z=S=\underset{\mathrm{0}}{\overset{z_t}{\int }}{\displaystyle \frac{\partial c}{\partial t}}\mathrm{d}z+\underset{\mathrm{0}}{\overset{z_t}{\int }}w{\displaystyle \frac{\partial c}{\partial z}}\mathrm{d}z.\end{array}$$

For our measurements, we choose z_t=125 m, so all biotic sources of carbon dioxide are within the column and S represents NEE.

Between a reference time t₀ (see below) and the time of a flight t_F the column has accumulated

$$\begin{array}{}\text{(9)}& {\displaystyle }\underset{t_{\mathrm{0}}}{\overset{t_{\mathrm{F}}}{\int }}S\mathrm{d}t& {\displaystyle }=\underset{t_{\mathrm{0}}}{\overset{t_{\mathrm{F}}}{\int }}\underset{\mathrm{0}}{\overset{z_t}{\int }}{\displaystyle \frac{\partial c}{\partial t}}\mathrm{d}z\mathrm{d}t+\underset{t_{\mathrm{0}}}{\overset{t_{\mathrm{F}}}{\int }}\underset{\mathrm{0}}{\overset{z_t}{\int }}w{\displaystyle \frac{\partial c}{\partial z}}\mathrm{d}z\mathrm{d}t,\text{(10)}& {\displaystyle }& {\displaystyle }=\underset{\mathrm{0}}{\overset{z_t}{\int }}c(z,t_{\mathrm{F}})\mathrm{d}z-\underset{\mathrm{0}}{\overset{z_t}{\int }}c(z,t_{\mathrm{0}})\mathrm{d}z+\underset{t_{\mathrm{0}}}{\overset{t_{\mathrm{F}}}{\int }}\underset{\mathrm{0}}{\overset{z_t}{\int }}w{\displaystyle \frac{\partial c}{\partial z}}\mathrm{d}z\mathrm{d}t.\end{array}$$

Dividing the accumulated amount of CO₂ by t_F−t₀ yields NEE averaged over this time span, denoted $\stackrel{\mathrm{‾}}{S}$:

$$\begin{array}{}\text{(11)}& \begin{array}{ccccc}\frac{{\displaystyle }{\int }_{t_{\mathrm{0}}}^{t_{\mathrm{F}}}S\mathrm{d}t}{{\displaystyle }{t}_{\mathrm{F}}-t_{\mathrm{0}}}& =& \stackrel{\mathrm{‾}}{S}& =& \frac{{\displaystyle }{\int }_{\mathrm{0}}^{z_t}c(z,t_{\mathrm{F}})\mathrm{d}z-{\displaystyle }{\int }_{\mathrm{0}}^{z_t}c(z,t_{\mathrm{0}})\mathrm{d}z}{{\displaystyle }{t}_{\mathrm{F}}-t_{\mathrm{0}}}\\ & & & & (\text{A)}\\ & & & +& \frac{{\displaystyle }{\int }_{t_{\mathrm{0}}}^{t_{\mathrm{F}}}{\int }_{\mathrm{0}}^{z_t}w{\displaystyle \frac{\partial c}{\partial z}}\mathrm{d}z\mathrm{d}t}{{\displaystyle }{t}_{\mathrm{F}}-t_{\mathrm{0}}}.\\ & & & & (\text{B)}\end{array}\end{array}$$

Term A represents the enhancement in CO₂ concentration and term B the vertical exchange of CO₂. We choose t₀ as the time when the surface radiation balance becomes negative, i.e. the time when the stable NBL starts to form. A positive $\stackrel{\mathrm{‾}}{S}$ means emission of CO₂ from the surface into the atmosphere.

The CO₂ concentration c can be calculated from CO₂ dry-air mole fraction $x_{{\mathrm{CO}}_{\mathrm{2}}}$, air temperature T and dry pressure p_d using the ideal gas constant R (Foken et al., 2012, p. 5):

$$\begin{array}{}\text{(12)}& c={\displaystyle \frac{p_{\mathrm{d}}}{RT}}{x}_{{\mathrm{CO}}_{\mathrm{2}}}.\end{array}$$

As COCAP measures $x_{{\mathrm{CO}}_{\mathrm{2}}}$, T, p and relative humidity U, the integral ${\int }_{\mathrm{0}}^{z_t}c(z,t_{\mathrm{F}})\mathrm{d}z$ in Eq. (11) term A is readily computed. However, each air parcel is sampled only once, at the time t_F when it passes the Fendt site. In order to evaluate the second integral in term A, ${\int }_{\mathrm{0}}^{z_t}c(z,t_{\mathrm{0}})\mathrm{d}z$, we assume horizontal and vertical homogeneity of the CO₂ dry-air mole fraction at the time t₀; i.e. $x_{{\mathrm{CO}}_{\mathrm{2}}}\left(t_{\mathrm{0}}\right)$ is assumed to be constant within the spatial domain relevant for our experiments. Thus we can calculate the second integral from the measurement of a different column at Fendt at t₀.

Note that this is a weaker assumption than the horizontal homogeneity of the CO₂ concentration in the NBL presumed in other studies (Choularton et al., 1995; Culf et al., 1999; Acevedo et al., 2004). All natural environments exhibit a certain horizontal heterogeneity in S. Daytime and nighttime CO₂ flux of a vegetated area are usually of the same order of magnitude, although different in sign. Before t₀, the convective boundary layer is well mixed up to a height of typically 1 km, whereas after t₀ a strong NBL confines emissions from the surface to the lowest ≈100 m of the atmosphere (see Sect. 4.3). Therefore, the horizontal heterogeneity in c caused by the horizontal heterogeneity in S is 1 order of magnitude smaller during the day than during the night. The convective mixing during the day also keeps vertical gradients inside the boundary layer low; hence the approximation of $x_{{\mathrm{CO}}_{\mathrm{2}}}\left(t_{\mathrm{0}}\right)$ being independent of z is justified.

Term B in Eq. (11) contains the product $w\partial c/\partial z$, which generally includes both turbulent exchange and subsidence. However, when a stable NBL has developed, little turbulent exchange takes place across the top of the NBL. In the statically neutral residual layer above the NBL, turbulence is present, but the vertical concentration gradient in the residual layer and as a consequence the net vertical transport of CO₂ is small. Hence, we neglect turbulent exchange and identify $w\partial c/\partial z$ with subsidence or lifting. The vertical wind speed w due to subsidence at a height of 100 m is usually on the order of 100 m d⁻¹, i.e. very low and therefore challenging to measure. We retrieve an estimate of w from the Integrated Forecast System (IFS) run by the European Centre for Medium-Range Weather Forecasts (ECMWF).

In order to calculate $\partial c/\partial z$ at different times between t₀ and t_F, we use a simple model for the growth of the NBL:

$$\begin{array}{}\text{(13)}& x_{{\mathrm{CO}}_{\mathrm{2}}}(t,z)=x_{{\mathrm{CO}}_{\mathrm{2}}}\left(t_{\mathrm{F}},{\displaystyle \frac{t_{\mathrm{F}}-t_{\mathrm{0}}}{t-t_{\mathrm{0}}}}z\right)(t_{\mathrm{0}}<t\le t_{\mathrm{F}}).\end{array}$$

Where $\frac{t_{\mathrm{F}}-t_{\mathrm{0}}}{t-t_{\mathrm{0}}}z$ exceeds the maximum height of the profile measured at time t_F, we assume the CO₂ dry-air mole fraction to be equal to $x_{{\mathrm{CO}}_{\mathrm{2}}}\left(t_{\mathrm{0}}\right)$. This model for the growth of the NBL can be visualised best by starting at t_F and looking back in time. At t=t_F the factor $(t_{\mathrm{F}}-t_{\mathrm{0}})/(t-t_{\mathrm{F}})$ is equal to unity and the model yields the measured profile. At earlier times, the measured profile is compressed in the z direction, such that the height of the NBL decreases linearly as we go back in time. As t approaches t₀, the model yields a thin layer enriched with CO₂ at the surface and a constant CO₂ dry-air mole fraction of $x_{{\mathrm{CO}}_{\mathrm{2}}}\left(t_{\mathrm{0}}\right)$ above.

The concentration c(t,z) is calculated from $x_{{\mathrm{CO}}_{\mathrm{2}}}(t,z)$ using Eq. (12). To this end, we determine p_d(t,z) and T(t,z) by linear interpolation in time between the first profile of the night and the profile measured at t_F.

In summary, the model for the growth of the NBL represents four simplifying assumptions: (1) during the night, the NBL height increases linearly; (2) the integral ${\int }_{\mathrm{0}}^{z_t}{x}_{{\mathrm{CO}}_{\mathrm{2}}}(t,z)\mathrm{d}z$ increases linearly with time; (3) the shape of the $x_{{\mathrm{CO}}_{\mathrm{2}}}$ profile within the NBL remains the same throughout the night; and (4) the dry pressure and temperature of an air column measured at Fendt are representative for the whole footprint of the measurement (see Sect. 3.4).

3.4 Footprint calculation

The columns of air probed at Fendt at different times had a different history, depending on the wind field and atmospheric stability. Atmospheric transport models can identify the surface areas that have contributed to an observed tracer concentration, i.e. the footprint of an observation. We simulate atmospheric transport with STILT, the Stochastic Time-Inverted Lagrangian Transport model (Lin et al., 2003; Gerbig et al., 2003a), which is based on NOAA's HYSPLIT particle dispersion model (Stein et al., 2015). In our configuration, STILT launches 10 000 air parcels at different heights (see below) at every full hour during the period of our NBL measurements. Driven by meteorological data with a resolution of 0.1^∘ × 0.1^∘ from the ECMWF IFS (European Centre for Medium-Range Weather Forecasts Integrated Forecast System), STILT calculates the back trajectories of these parcels until 10 h in the past. For each time step of the simulated transport, the model determines the sensitivity of the CO₂ concentration in the parcel to the CO₂ flux at the surface. To do so, the height up to which mixing occurs is estimated from the meteorological data using a modified Richardson number method (Lin et al., 2003). Surface fluxes influence air parcels within a column that extends from the surface to one-half this height in each time step (Gerbig et al., 2003b).

The back trajectories calculated by STILT are then aggregated into mole fraction footprints on a regular grid with a resolution of 2 km×2 km. As explained in the previous section, we assume the $x_{{\mathrm{CO}}_{\mathrm{2}}}$ distribution to be homogeneous in the lateral and horizontal directions at time t₀. We therefore restrict the aggregation to that part of each back trajectory that lies between t₀ and the time of measurement.

A single STILT run determines the sensitivity of an observation at a specific height to upwind fluxes. Formally, the mole fraction footprint of a measurement taken at the geographic location (l₁,l₂) at time t and observation height z can be written as $f(l_{\mathrm{1}},l_{\mathrm{2}},z,t\mid l_{G\mathrm{1}},l_{G\mathrm{2}})$. As all our measurements were taken at the same horizontal location, the dependency of f on l₁ and l₂ will be omitted hereafter. The mole fraction footprint is a function whose value is the sensitivity to the surface flux at the grid cell specified by (l_G1,l_G2) in units of $\left[f\right]=\mathrm{\mu }\mathrm{mol}{\mathrm{mol}}^{-\mathrm{1}}\mathrm{\mu }{\mathrm{mol}}^{-\mathrm{1}}{\mathrm{m}}^{\mathrm{2}}\mathrm{s}$. To determine the relative contribution of surface fluxes in different areas to our NBL-derived fluxes, we need a different but related function, the flux footprint f_F with units $\left[f_{\mathrm{F}}\right]=\mathrm{\mu }\mathrm{mol}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{s}}^{-\mathrm{1}}\mathrm{\mu }{\mathrm{mol}}^{-\mathrm{1}}{\mathrm{m}}^{\mathrm{2}}\mathrm{s}=\mathrm{1}$. The flux footprint is calculated by integration over an array of mole fraction footprints for different measurement heights, i.e. analogous to Eq. (11) term A and Eq. (12):

$$\begin{array}{}\text{(14)}& f_{\mathrm{F}}\left(t\right|l_{\mathrm{G}\mathrm{1}},l_{\mathrm{G}\mathrm{2}})={\displaystyle \frac{{\displaystyle }\underset{\mathrm{0}}{\overset{z_t}{\int }}{\displaystyle \frac{p_{\mathrm{d}}(z,t)}{RT(z,t)}}\cdot f(z,t|l_{\mathrm{G}\mathrm{1}},l_{\mathrm{G}\mathrm{2}})\mathrm{d}z}{t-t_{\mathrm{0}}}}.\end{array}$$

Dry pressure p_d and air temperature T at time t and height z are inter- or extrapolated from the measured profiles. The ensemble of mole fraction footprints comprises footprints for 12 different measurement heights between 10 and 120 m in 10 m steps.

The meteorological data we use have a horizontal resolution of 0.1^∘ × 0.1^∘, corresponding to 11 km×8 km at the latitude of Fendt. Terrain features that are smaller than a grid cell, like the valley slope to the west of the Fendt site, cannot be represented at these resolutions. The vertical resolution of the meteorological data depends on height above ground. The lowest layer extends from the ground to 10 m height, and the following five layers extend from the top of the previous layer to 31, 55, 80, 108 and 138 m. The temporal resolution of the ECMWF IFS data is 3 h.

4 Results and discussion

4.1 Uncertainty of $x_{{\mathrm{CO}}_{\mathrm{2}}}$ measurements

The uncertainty of COCAP's $x_{{\mathrm{CO}}_{\mathrm{2}}}$ measurements due to drift and calibration errors is about 1 µmol mol⁻¹ (Kunz et al., 2018). The additional uncertainty caused by noise is dependent on the data treatment, as can be seen from Fig. 5. This Allan deviation plot (Allan, 1987) is based on measurements of a gas standard ($x_{{\mathrm{CO}}_{\mathrm{2}}}=\mathrm{447.44}{\text{µmol}\cdot \text{mol}}^{-\mathrm{1}}$) over a period of 1.4 h, taken in the field on 6 July 2016.

Figure 5 Allan deviation of CO₂ dry-air mole fraction σ_xCO2 as measured after deconvolution and after deconvolution followed by smoothing (convolution with a Gaussian function of 5 s FWHM) for different averaging periods τ. All three cases converge for averaging periods longer than 100 s.

The curves illustrate that deconvolution amplifies noise in the data by a factor of 7 if no averaging is applied (τ=1 s). However, if more than 100 samples are averaged (τ≥100 s), the difference between original and deconvoluted data becomes negligible and the uncertainty of the average due to noise is lower than 0.5 µmol mol⁻¹. All our column integrals (see Sect. 3.3) have a sample size larger than 100.

Figure 5 also shows that the Allan deviation of deconvoluted data that have been smoothed by convolution with a Gaussian function of 10 s full width at half maximum (FWHM) increases between τ=1 s and τ=5 s. This increase is an artefact caused by the autocorrelation that the smoothing induces. If COCAP were perfectly calibrated and exhibited no drift, any single point in the smoothed data set would have an uncertainty of 2.1 µmol mol⁻¹ (corresponding to τ=5 s), not 0.8 µmol mol⁻¹ (corresponding to τ=1 s).

4.2 Meteorological conditions

From the data collected during ScaleX 2016 we calculate NEE for the nights 6–7 July and 9–10 July. The sun set at 19:15 and 19:14 UTC on 6 and 9 July, respectively, and rose at 03:25 and 03:27 UTC on 7 and 10 July, respectively. Both nights were free of precipitation. Cloud cover was high during the first night (see Fig. 6a), but the pronounced negative net radiation (Fig. 6b) indicates that the clouds were mostly transparent for outgoing long-wave radiation.

Figure 6 Meteorological conditions during the NBL soundings: (a) cloud cover, (b) net radiation E^*, (c) temperature T of air (2 m height) and soil (2 cm depth), and (d) horizontal wind speed u (3.5 m height). Cloud cover was determined at MOHP, and all other observations were made directly at Fendt. Time is given in Coordinated Universal Time.

On the second night the sky was clearer, resulting in a steadier radiation balance. During both nights, strong radiative cooling was observed. Air temperature decreased from 18 to 9 ^∘C and from 24 to 11 ^∘C over the course of the first and second nights, respectively (Fig. 6c). In combination with low wind speeds (Fig. 6d) this led to the development of a pronounced temperature inversion at the surface, i.e. a stable NBL. The change from positive to negative net radiation occurs approximately at t₀=18:00 UTC on both nights.

4.3 Profiles

We carried out a total of 27 flights during the ScaleX 2016 campaign. For the calculation of a NBL budget, we analyse those flights that took place after t₀=18:00 UTC and reached a height of at least 125 m. Twelve flights fulfil these criteria: flights 4 through 10 (first night, Fig. 7) and flights 19 through 23 (second night, Fig. 8). For display in panel (b) of Fig. 7 and 8 the CO₂ dry-air mole fraction has been smoothed with a Gaussian filter of 5 s FWHM. To prevent distortion in the vertical direction, the height above ground level z has been filtered the same way. For this reason, the upper end of the profiles in panel (b) is at a slightly lower height than in panel (a). Calculation of the NBL fluxes (see Sect. 4.5) was carried out with unfiltered $x_{{\mathrm{CO}}_{\mathrm{2}}}$ and z.

The times given in Figs. 7 and 8 are the midtimes of the flights rounded to a full 10 min for readability. The exact times of take-off and landing are provided in Table 2.

Table 2 Take-off and landing times (UTC).

Figure 7 Profiles from flights 4 through 10, carried out in the night from 6 to 7 July 2016. Times are given in Coordinated Universal Time and specify the middle of the flight, rounded to the next 10 min. (a) Virtual potential temperature θ_v at height above ground level z. (b) Carbon dioxide dry-air mole fraction $x_{{\mathrm{CO}}_{\mathrm{2}}}$ at height above ground level z. To reduce noise $x_{{\mathrm{CO}}_{\mathrm{2}}}$ has been smoothed (for details see text). In both (a) and (b) the light grey curves are copies of the previous profile (ascent, horizontal flight and descent combined). (c) Flight track with horizontal projections (grey) for clearness. The northward and eastward directions are marked. The tick marks at the ground plane are 100 m apart. The location of take-off and landing differs slightly between flights, but all flights took place within the same 250 m × 250 m × 150 m bounding box.

Figure 8 Same as Fig. 7, but for flights 19 through 23 carried out in the night from 9 to 10 July 2016. The axis for virtual potential temperature is shifted towards higher temperatures compared to Fig. 7 but covers the same span. All other axes are unchanged. The horizontal legs were flown at 10 m height.

During the first night, a stable NBL can be identified from the UAS profiles for flights 6 through 10. The upper end of the temperature inversion aligns with the top of the CO₂ enhancement to within 10–20 m. At the time of flights 6 and 8 through 10, the NBL has a height of 50–70 m, whereas the profile from flight 7 indicates a greater NBL height of ≈100 m. We interpret this as an indication that the column measured in flight 7 has been influenced by katabatic inflow of cool, CO₂-enriched air at some point during the night, potentially hours before the flight and kilometres away from Fendt. This interpretation is supported by the flux estimates (see Sect. 4.5). The profiles from flight 5, which exhibit virtually no gradient, are discussed below.

During the second night, a stable NBL with a height of 50–70 m is visible in all profiles. The flight pattern had been refined and included two ascents and descents far enough from each other to avoid disturbance of the measurements in the second part by air movements caused during the first part. These redundant measurements give insight into the reliability of the measurement system and the variability of temperature and CO₂ dry-air mole fraction on small temporal and spatial scales. The data from flight 21 agree well between each of the two ascents and descents, suggesting that disturbances by the UAS, instrument noise and drift are small compared to the observed signals. Flights 22 and 23 were carried out only 1 and 2 h later, respectively, and followed the same flight track. However, the data from these flights reveal considerable differences between each of the two ascents and descents, especially in $x_{{\mathrm{CO}}_{\mathrm{2}}}$ for heights below 50 m. We interpret this as natural variability on the scale of the flight track, i.e. ≈200 m in horizontal distance and ≈3 min in time. This small-scale variability is a source of random error in NBL budgets. In our flux calculations multiple ascents during the same flight are effectively averaged, resulting in a reduction of the random error.

The $x_{{\mathrm{CO}}_{\mathrm{2}}}$ profiles measured during flights 20, 22 and 23 all exhibit a non-zero gradient with height in the region above the strong inversion, indicating that some CO₂ has escaped the stable NBL. This is supported by the profiles of virtual potential temperature, which are more inclined above the NBL in comparison to the first night. Both features might be the result of intermittent turbulence, a phenomenon often observed at night that can have different causes (see Aubinet, 2008, and references therein). Our budgets include the measurements up to 125 m height, so any CO₂ that has been transported higher than this is missing in the budgets. In future campaigns, flights with a greater maximum height could be carried out to quantify the effect this has on the NBL-derived fluxes, or to extend the budget vertically.

The flight pattern used during the second night also included two horizontal transects at 10 m height that were flown at a ground speed of 3 m s⁻¹. Their purpose was to enable measurements of undisturbed air near the ground, but later analysis of flight 14 (see Sect. 3.1) revealed that the ground speed was insufficient to fully reach this goal.

The profiles for flight 5 are close to straight vertical lines, which would indicate a well-mixed atmosphere. However, they were measured under low wind speed 1 full hour after the surface radiation balance became negative, i.e. under conditions favourable for the development of a stable nocturnal boundary layer and accumulation of CO₂ near the ground. This apparent contradiction can be explained by comparing COCAP's data to tower-based measurements. Figure 9a shows the CO₂ profile taken by COCAP together with data from the 9 m mast and from HPB (see Sect. 2.1 and 2.2).

Figure 9 Comparison of $x_{{\mathrm{CO}}_{\mathrm{2}}}$ measurements by COCAP, by the ICOS station HPB (a point measurement 460 m above ground level at Fendt as described in Sect. 2.2, representing the CO₂ dry-air mole fraction in the residual layer) and by the on-site 9 m mast (a) for flight 5 and (b) for flight 8. The dots and bars are the mean and standard deviation, respectively, for each 1 min sampling period of the mast within the time interval from 15 min before take-off to 15 min after landing. Note that the scaling of the vertical axis changes at height z=10 m.

The diagram includes those measurements from the mast that fall into the time interval from 15 min before take-off to 15 min after landing. They reveal that the CO₂ dry-air mole fraction near the ground was increased relative to the upper two-thirds of the profile and fluctuated strongly, e.g. between 450 and 650 µmol mol⁻¹, at 3 m height. These observations are in line with a weakly stable layer near the surface: surface fluxes accumulated in this layer, but weak turbulent events caused by wind shear, for example, occasionally spread them out to higher layers. The disturbance by the multicopter during take-off or landing prevented COCAP from capturing this accumulation. On the other hand, the higher part of COCAP's profile, taken in the residual layer that is left over from the daytime mixed layer, matches the mean CO₂ dry-air mole fraction measured at HPB during the time interval from 1 h before take-off to 1 h after landing. This agreement confirms that COCAP was working properly during the flight.

The profile from flight 8, carried out later in the same night, is consistent with the measurements at the 9 m mast (Fig. 9). We see two reasons for this difference to flight 5. Firstly, the radiative cooling (see Fig. 6) at the time of flight 8 (23:10) was stronger than at the time of flight 5 (19:10). The temperature gradient near the ground was not resolved during flight 5, but the weaker radiative cooling compared to the later flight has likely resulted in a weaker temperature inversion that allowed more vertical displacement of air by the multicopter. Secondly, the thicker NBL at 23:10 with a less steep CO₂ gradient close to the ground means that potential sampling of air parcels originating from above or below the multicopter did not affect the measurements during flight 8 as much as during flight 5.

At heights above 70 m the CO₂ profile from flight 8 approaches the measurements at HPB, indicating that the stable NBL retains most of the surface-emitted CO₂. Likewise, the CO₂ profiles of all other flights come near the measurements at HPB above the NBL (not shown). This suggests that any transport of CO₂ across the top of the NBL is small in magnitude.

The measurement of continuous profiles of the CO₂ dry-air mole fraction up to heights of 100 m or more has been challenging in the past. In some studies, NBL budgets were therefore based on a measurement near the ground and an assumed gradient up to the top of the NBL. However, the complex shape of the profiles displayed in Figs. 7 and 8 suggests that neither the assumption of a constant (see Acevedo et al., 2004) nor a linearly decreasing (see Culf et al., 1997) CO₂ dry-air mole fraction would properly represent the conditions at Fendt. The detailed structures resolved in our measurements also indicate great potential of combined measurements of meteorological parameters and trace gas mole fractions for studying small-scale phenomena in the NBL.

4.4 Disturbance by the UAS

The potential virtual temperature measured at heights between 10 and 60 m is generally higher during descent than during ascent. This effect is more pronounced for flights 19 through 23 (Fig. 8), likely due to the stronger temperature gradient compared to flights 5 through 10 (Fig. 7). The observed difference supports the reasoning of Sect. 3.1: As the multicopter descends, the onboard sensors measure warmer air that was pushed downwards by the rotors. Close to the ground (at heights below 10 m) closed flow loops start to form and colder air from below the multicopter reaches the sensors during descent, as can be seen in the profiles from flights 6, 8, 9, 10 and 23.

Systematic differences between ascent and descent are less visible in the profiles of CO₂ dry-air mole fraction, likely due to a larger variability of CO₂ within the nocturnal boundary layer. This variability is reflected in the difference in $x_{{\mathrm{CO}}_{\mathrm{2}}}$ between each of the two ascents and descents in the flights 19 through 23, especially flights 20 and 23.

Flight 14 was dedicated to the investigation of vertical mixing during horizontal movement at different air speeds. It was carried out on 7 July at 22:15 UTC. Winds were particularly low that night (on average 0.3 m s⁻¹ between 22:00 and 22:30 UTC), and hence ground speed of the UAS was approximately equal to air speed. A stable nocturnal boundary layer had developed, as can be seen from the profiles of θ_v and $x_{{\mathrm{CO}}_{\mathrm{2}}}$ measured during an earlier flight at 20:15 UTC (see Fig. 10a and b).

Figure 10 Vertical profiles of (a) virtual potential temperature θ_v and (b) CO₂ dry-air mole fraction $x_{{\mathrm{CO}}_{\mathrm{2}}}$ at 20:15 UTC on 7 July 2016. (c) Track of flight at 22:15 UTC on the same night coloured by horizontal ground speed. Three sections of nominal speeds of 5, 3 and 1 m s⁻¹ are marked with brackets. Height was 10 m above ground throughout the flight. (d) Median of virtual potential temperature and CO₂ dry-air mole fraction measured during those three sections of the flight. Bars represent the bootstrapped standard error of the median; see text for details. The standard error of the virtual potential temperature is so small that the vertical bars are barely visible. At lower speeds, θ_v is lower and $x_{{\mathrm{CO}}_{\mathrm{2}}}$ is higher, suggesting the sampling of air that originates from below the flight height.

The UAS flew a spiral pattern at a height of 10 m above ground with decreasing ground speed (Fig. 10c). Throughout flight 14, COCAP's air inlet faced the direction of movement. The flight took place over a flat, homogeneous meadow. Hence, we assume that terrain and vegetation had caused no heterogeneity in the lateral distribution of temperature and CO₂. We analyse three sections of nominal speeds of 5, 3 and 1 m s⁻¹. Figure 10d shows the median virtual potential temperature and CO₂ dry-air mole fraction for each section. The standard error of the median was calculated by bootstrapping with 1000 samples generated from the empirical distribution of the measurements (Wilcox, 2012, pp. 43) and is depicted as horizontal and vertical bars.

The decrease in virtual potential temperature with decreasing speed in Fig. 10d suggests that upward mixing of air from lower layers has a stronger influence on the measurements at lower speed. Likewise, the CO₂ dry-air mole fraction measured at 1 m s⁻¹ is 20 µmol mol⁻¹ higher than during faster flight. However, we did not observe a significant difference in $x_{{\mathrm{CO}}_{\mathrm{2}}}$ between a ground speed of 3 and 5 m s⁻¹. The sample inlet for the CO₂ measurement extends 20 cm to the side of the rotors, while temperature and humidity are measured directly below a rotor (see Fig. 3). As the sample inlet was pointing forward throughout the flight, it might have mostly avoided partially closed flow loops during movement at 5 m s⁻¹, while the temperature and humidity sensors were still affected.

In summary, our results suggest that measurements taken during the ascent of the multicopter are more reliable than those taken during descent and hover. Horizontal transects at low heights can yield measurements that are contaminated with air from below the sampling height. This contamination is lower at higher horizontal air speed, because the multicopter moves away from the vortices it has created. Our experiment does not answer the question of whether a horizontal speed of 5 m s⁻¹ at 10 m height is sufficient to avoid the contamination entirely.

4.5 Carbon dioxide fluxes

The first profiles of the first and second nights were taken at 18:10 UTC (flight 4) and 20:10 UTC (flight 19), respectively. Hence, flight 4 is representative for the $x_{{\mathrm{CO}}_{\mathrm{2}}}$ profile at $t_{\mathrm{0}}=\mathrm{18}:\mathrm{00}$ UTC, but flight 19 is not. We therefore need an estimate for the profile at t₀. Due to the convective mixing that takes place during the day, the CO₂ dry-air mole fraction within the boundary layer is nearly independent of height, an assumption that is supported by the profile from flight 4 (see Fig. 7). Assuming further that all surface fluxes were trapped in the developing NBL, air parcels above the NBL height should have preserved the CO₂ dry-air mole fraction of the column between t₀ and the time of the first flight. Consequently, we assume the whole column $x_{{\mathrm{CO}}_{\mathrm{2}}}(t_{\mathrm{0}},z)$ to be equal to the mean dry-air mole fraction of the first measured profile between 50 and 125 m height. For consistency we apply this approach to both nights.

The fluxes we calculated from the NBL budgets are listed in Table 3, given as the amount of CO₂ per time and surface area.

Table 3 Fluxes of CO₂ calculated from NBL budgets. Begin and end times are given in UTC. The end time is specified as the midtime of the portion of the flight used for determination of the NBL budget.

The storage flux in Table 3 corresponds to term (A) in Eq. (11), the subsidence flux corresponds to term (B) and the total flux is equal to $\stackrel{\mathrm{‾}}{S}$, i.e. the NEE averaged over the time from t₀ to t_F. During both nights, horizontal convergence of air masses led to lifting and consequently a negative subsidence flux. However, the subsidence flux was small compared to the storage flux, accounting for about 1 % of the total flux. An important consequence of the low subsidence flux is that errors stemming from the simplified model of the NBL growth (see Sect. 3.3) have only a minor influence on the uncertainty of the total flux.

The plausibility of our results can be checked against EC and chamber measurements taken at Fendt. Both the EC and the chamber measurements observed only the fluxes from the pasture at the site, while the NBL budget has a larger footprint. Even at low wind speeds of 0.5 m s⁻¹ air parcels travel 1.8 km every hour. Therefore, the NBL budget also includes sources that are located several kilometres apart. Given the land cover around Fendt, those sources likely include forests, crop fields and potentially some residential areas (see Figs. 14 and 15 for exemplary footprints). Nevertheless, as pasture is the dominant land cover in the area, all three methods should agree on the order of magnitude of the CO₂ flux at night.

No EC measurements of acceptable quality are available for either of the nights we probed the NBL (Fig. 11). Conditions of strong radiative cooling combined with weak wind resulted in stable conditions and a violation of the assumptions underlying the EC technique. As a backup, we calculated the mean diurnal cycle from the EC measurements taken between 4 July 2016, 00:00 UTC, and 11 July 2016, 23:59 UTC, a period that includes all our flights and was reasonably consistent in the diurnal variations in temperature. The result is presented in Fig. 11. All fluxes calculated from the NBL budget lie within the range of NEE observed by EC between 18:00 and 01:00 UTC (6–16 $\mathrm{\mu }\mathrm{mol}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{s}}^{-\mathrm{1}}$). The later the flight at Fendt took place, the lower the NBL-based average NEE, indicating a decreasing flux over the course of the night. We interpret this, at least partially, as an effect of the temperature decrease during the night (Fig. 6), which reduces respiration. In contrast, NEE measured by the EC station increases during the night. However, an increase in respiration over the course of the night is implausible. Given the small number of EC measurements of acceptable quality, this apparent trend is likely an artefact.

Figure 12 shows the NBL-derived fluxes in comparison to chamber measurements. Data from the small chambers are available only for the second night. Opaque chambers measure respiration, while clear chambers, the EC station and the NBL budget observe NEE. Therefore, a comparison of the fluxes obtained with these different techniques is only meaningful when photosynthesis is low or absent, i.e. roughly between sunset and sunrise. The convergence of the fluxes of the clear and dark chambers just after 18:00 UTC suggests that photosynthesis has largely ceased as early as t₀. Hence, throughout the time span for which we determine the NBL budget NEE is dominated by respiration, and all the different techniques are comparable. Surprisingly, the measurements with the big chambers yield fluxes only one-third as high as obtained with the small chambers, even though all chambers were deployed close to each other on the same meadow. Despite careful investigation, the reason for this discrepancy has not yet been found. The NBL budget agrees in magnitude to the fluxes measured with the small chambers. Similarly to the NBL budget, all chamber measurements exhibit a negative trend in fluxes over the course of both nights.

Figure 11 Comparison of vertical CO₂ fluxes F calculated from the NBL budget and using the EC method. The fluxes from the NBL budget are depicted as lines, where the vertical position of each line represents the average flux over the time span specified by the horizontal extent of the line. Open circles represent the quality-filtered EC measurements taken on the same days as the NBL measurements. Solid dots represent the mean diurnal cycle of the quality-filtered EC measurements averaged over the period from 4 July 2016, 00:00 UTC, to 11 July 2016, 23:59 UTC. Upward and downward arrows mark the time of sunrise and sunset, respectively. NBL and EC agree in magnitude of NEE at night but not in sign of trend.

Figure 12 Comparison of CO₂ flux F calculated from the NBL budget and flux measured with different chambers. The axis scaling and sunset–sunrise markers are the same as in Fig. 3. Nighttime fluxes observed with the small and the big chambers differ by a factor of 3 for unknown reasons. The NBL budget agrees in magnitude and sign of trend to the measurements with the small chambers.

In addition to in situ measurements at Fendt, the range of nighttime NEE of pasture and forests observed in other studies at central European sites with a climate similar to Fendt (Cfb or Dfb in the Köppen–Geiger classification according to Peel et al., 2007) provides a plausibility check for the NBL budgets (Table 4).

Gilmanov et al. (2007)Wohlfahrt et al. (2005)Wohlfahrt et al. (2005)Chen et al. (2014)Siebicke et al. (2012)Longdoz et al. (2008)Kutsch et al. (2008)Kutsch et al. (2008)Hommeltenberg et al. (2014)

Table 4 Nighttime CO₂ fluxes observed in other studies (minimum–maximum of reported values).

We exclude crop fields from the comparison, as their NEE depends heavily on crop type and time of harvest. Compared to the literature values, NEE for Fendt derived from the NBL budget is on the high end of ranges reported for pasture and higher than most fluxes reported for forests. One explanation is that our measurements took place on two fair weather days in the warmest month of the year 2016, which likely resulted in higher respiration than observed on average over a longer period. Furthermore, Fendt lies in a region with organically rich soils (Fig. 2a). Soil organic carbon content has been shown to be positively correlated with microbial biomass (Habashi, 2016), suggesting particularly strong respiration under beneficial conditions. This explanation is supported by the measurements at Mooseurach (Table 4), a drained peatland forest 20 km to the east of Fendt, where respiration fluxes of up to 15 $\mathrm{\mu }\mathrm{mol}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{s}}^{-\mathrm{1}}$ have been observed.

Another potential cause for higher fluxes observed with the NBL budget relates to the terrain at Fendt. At night, katabatic flows of cool, CO₂-rich air can stream down the steep slope west of the measurement site. Though Fendt is situated in a valley with only a shallow slope to the east, this inflow might lead to localised lifting of air that is not accounted for in the ECMWF IFS data and hence not included in our calculation of subsidence. The increased NBL height and high variability in the lowest 50 m observed during flight 7 as well as the higher flux derived from the NBL budget are an indication of such an inflow event.

4.6 Sensitivity of fluxes

The NBL budget is influenced by measurement uncertainty, incomplete knowledge about the state of the atmosphere and data selection. In order to quantitatively assess the influence of these factors on our results, we changed the procedure of calculating fluxes in one of the following ways:

1. by adding a bias of ±2 m to the altitude measurements;

2. by adding a bias of $\pm \mathrm{3}\mathrm{\mu }\mathrm{mol}{\mathrm{mol}}^{-\mathrm{1}}$ to $x_{{\mathrm{CO}}_{\mathrm{2}}}$ of all but the first profile of each night;

3. by using COCAP data for the whole column instead of replacing $x_{{\mathrm{CO}}_{\mathrm{2}}}$ in the lowest 9 m with measurements taken at the 9 m mast;

4. by using COCAP data taken during the whole flight, i.e. using ascent, descent and hover instead of ascent only; or

5. by disregarding subsidence.

Check 1 accounts for the uncertainty of COCAP's pressure-based altitude measurements. Check 2 allows us to evaluate the influence of both the uncertainty of COCAP's $x_{{\mathrm{CO}}_{\mathrm{2}}}$ measurements and the spatial heterogeneity of x(t₀). The former is known from experiment (see Sect. 4.1), and the latter can be estimated from the CO₂ measurements at HPB. Assuming that the 131 m inlet at HPB is in the residual layer all night, the interquartile range of the $x_{{\mathrm{CO}}_{\mathrm{2}}}$ measurements of a single night approximately reflects the variability of the background onto which fluxes accumulate. The interquartile range amounts to 1.1 and 2.4  µmol mol⁻¹ for the period from 18:00 to 02:00 UTC on the first and second nights, respectively. Checks 3 and 4 relate to the disturbance caused by the UAS, which is discussed in Sects. 3.1 and 4.4.

The mean fluxes for each night obtained using the changed procedures are summarised in Table 5.

Table 5 Sensitivity of CO₂ flux to different factors (see text for details). Fluxes are in micromoles per square metre per second.

The largest difference to the normal (“no change”) procedure occurs when $x_{{\mathrm{CO}}_{\mathrm{2}}}$ is altered (±11 % for the first and ±7 % for the second night). Changing the altitude or not using the data from the 9 m mast also has a considerable influence on the mean flux.

Figure 13 shows the values from Table 5 in graphical form. In addition, the fluxes calculated for each flight are depicted, visualising how their spread is affected by the different checks.

Figure 13 Sensitivity of CO₂ fluxes calculated from the NBL budgets to changes in calculation procedure. Dots denote the fluxes calculated for each flight, and crosses mark the mean of these fluxes. The lines are visual aids to facilitate comparison to the mean flux of the normal (“no change”) procedure.

A substantial increase in spread is observed only when the data from the 9 m mast are not used.

Overall, the results from the sensitivity checks indicate that the NBL method is robust against measurement uncertainty in the altitude and $x_{{\mathrm{CO}}_{\mathrm{2}}}$ measurements, spatial heterogeneity of x(t₀), disturbance of the NBL caused by the UAS, and the effect of subsidence. It should be noted that the mean vertical wind extracted from the ECMWF IFS model was relatively small during the two nights of our measurements. Under different conditions, e.g. in a strong high-pressure system, the effect of subsidence or lifting on the NBL budget could be much higher.

4.7 Flux footprint

Example flux footprints of one NBL budget of each night are visualised in Figs. 14 and 15.

Figure 14 Footprint of an NBL budget at Fendt on 6 July, 21:00 UTC. (a) Relative contribution of each grid cell to the total sensitivity of the budget to surface fluxes. Water bodies depicted for orientation. (b) Contour of all grid cells with a relative sensitivity of 1 % or higher on top of a simplified land cover map (CORINE 2012 v18.5, European Environment Agency, EEA (2016)). The area observed is dominated by forests, pasture and cropland.

Figure 15 Same as Fig. 14, but for 9 July, 21:00 UTC.

The footprint depicted in Fig. 14 was calculated for a column of air passing Fendt on 6 July 2016 at 21:00 UTC, i.e. close to the time of flight 6. The 1 % contour of the footprint encloses an area of 60 km², which accounts for 60 % of the total sensitivity. The land cover map suggests that the NBL budget represents mainly the respiration of forests, pasture and croplands north of Fendt, with little contribution from urban areas.

The footprint depicted in Fig. 15 was calculated for a column of air passing Fendt on 9 July 2016 at 21:00 UTC, i.e. close to the time of flight 20. The 1 % contour of the footprint encloses an area of 80 km², which accounts for 70 % of the total sensitivity. Again, the NBL budget is mainly influenced by forests, pasture and croplands.

The footprints for other times during the two nights are similar in size, i.e. on the order of 100 km². They mostly cover the sector within 20 km north-west to north-east of Fendt.

We recognise that the relatively low spatial and temporal resolution of the ECMWF IFS meteorological model entails errors in the transport modelling. Variability of the horizontal wind component within a grid cell and on timescales below 3 h is neglected, possibly resulting in an underestimation of the footprint size. Likewise, terrain features that are smaller than a grid cell are not represented in the meteorological model. However, as our NBL budgets cover timescales of 1–7 h and the footprints extend over many grid cells, sub-scale variability should play only a minor role. We are therefore confident that the model results provide a reasonable estimate of the region seen by the NBL budget method.

5 Conclusions and outlook

To the best of our knowledge, we have for the first time determined nocturnal boundary layer budgets based on trace gas measurements with an unmanned aircraft. During two nights we repeatedly sampled the NBL with a multicopter carrying COCAP, a lightweight analyser designed for deployment on unmanned aircraft. Simultaneous measurement of CO₂ dry-air mole fraction, air temperature, humidity and pressure allowed the quantification of the rate of accumulation of carbon dioxide in the NBL. By applying deconvolution we could improve the temporal resolution of the CO₂ measurements, thus achieving a higher vertical resolution of the profiles. We estimated the effect of subsidence or lifting on the NBL budgets with the help of weather forecast data and corrected the budgets accordingly. The respiration fluxes obtained from the NBL budgets are plausible in comparison to other flux measurements at the Fendt site, though on the high end of the range reported in the literature for sites with land cover and climatic conditions similar to Fendt. A potential positive bias in the obtained fluxes could be caused by convergence of cool, CO₂-rich air at the floor of the valley in which Fendt is located. The current data set does not allow us to confirm or rule out this effect. In a future campaign, however, simultaneous deployment of a second UAS on the elevated plateau west of the site could provide more insight, as downward transport of CO₂ should result in consistently lower accumulation and hence lower flux estimates obtained from NBL budgets on the plateau.

We have investigated how the disturbance of the NBL caused by a multicopter influences in situ measurements. We found that while flying close to the ground, air from below the UAS can reach the sensors, causing a bias if the respective quantity has a non-zero gradient. To prevent this bias from affecting the NBL budget we replaced the airborne $x_{{\mathrm{CO}}_{\mathrm{2}}}$ measurements taken at low height with measurements from a 9 m mast. At greater height, some of our profiles exhibit a systematic difference between ascent and descent. During descent, the airborne sensors are moved into a volume of air that may have been disturbed by the downwash of the multicopter's rotors. Therefore, we use only data captured during ascent for NBL budgeting.

The robustness of our approach has been demonstrated by a sensitivity analysis. The largest uncertainty of the NBL budget is caused by spatial heterogeneity of the CO₂ dry-air mole fraction in the late afternoon combined with the uncertainty of the CO₂ measurement. The estimated combined error in $x_{{\mathrm{CO}}_{\mathrm{2}}}$ results in ±11 % change of the mean of the fluxes obtained from the NBL budget for the first night. Using only data from the UAS and not from the 9 m mast increased the spread of the fluxes but changed their mean by no more than 10 %. This suggests that satisfactory NBL budgets can be determined from UAS data even if no stationary measurements near the ground are available. For future studies, we suggest positioning the sample inlet 50–100 cm above the rotors to further reduce the sampling of air that was displaced or mixed by the UAS.

The region that influences the NBL budget has often not been reported in past studies. We improved on this situation by carrying out mesoscale modelling. While the driving meteorological data and the underlying topography do not resolve small structures at and below the scale of 1 km, our method gives at least an estimate of the region that influences the NBL budget. Under the conditions of our measurements the footprints were on the order of 100 km² in size. In situ wind measurements would enable validation of the meteorological data and possibly improvement of the transport modelling. Such measurements could be taken by UASs without the need for additional sensors (Mayer et al., 2012; Neumann and Bartholmai, 2015).

Future NBL studies could employ multiple UASs simultaneously to quantify spatial heterogeneity and horizontal gradients in the CO₂ dry-air mole fraction. Firstly, this would support the analysis of the uncertainty of the NBL-derived fluxes. Secondly, concurrent profiles could yield constraints for the net advection of CO₂.

While we carried out our measurements with multicopters, fixed-wing aircraft would also be capable platforms for NBL studies. The vortices generated by their wings are slower and spread out wider than the concentrated downwash produced by the rotors of a multicopter. Therefore they should cause less interference with the NBL soundings and could provide precise measurements down to ground level. Additionally, their typically higher horizontal speed makes it easier to evade any disturbance that they create.

Another possibility to reduce the disturbance of measurements near the ground would be a different placement of the inlet. Given the asymmetric flow pattern below and above a multicopter's rotors (see Sect. 3.1), sampling from several rotor diameters above the UAS should reduce the artefacts caused by closed flow loops.

NBL budgets based on UAS measurements are an effective and efficient tool for the quantification of nocturnal fluxes. Besides ecosystem respiration, this tool could also be applied to detect carbon dioxide emissions of other sources, e.g. urban areas. Small and lightweight sensors for other tracers such as methane would open up even more possibilities. Alternatively, compact time-resolved sampling systems (Andersen et al., 2018) or long flexible tubing (Brosy et al., 2017) can be used in connection with conventional ground-based instrumentation to measure a whole range of species.

In summary, we have demonstrated that nocturnal surface flux estimates can be derived from UAS-based gas measurements by means of an NBL budget approach. Given the moderate cost of UASs and their minimal infrastructure requirements, this innovation makes the NBL budget method for the quantification of surface fluxes much more accessible. Spurred by the increasing adoption of unmanned aircraft in geoscience and the development of miniaturised high-accuracy sensors for different tracers, we foresee wide adoption of this technique in the coming years.