# Methodological development on conditional sampling method (part I): sensibility to statistical and technical characteristics

A. Fotiadi, Fabienne Lohou, A. Druilhet, Dominique Serça, Y. Brunet, Robert Delmas

► To cite this version:

A. Fotiadi, Fabienne Lohou, A. Druilhet, Dominique Serça, Y. Brunet, et al.. Methodological development on conditional sampling method (part I): sensibility to statistical and technical characteristics. Boundary-Layer Meteorology, 2005, 114, pp.615-640. 10.1007/s10546-004-1080-9 . hal-00138884

HAL Id: hal-00138884

https://hal.science/hal-00138884v1

Submitted on 27 May 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
METHODOLOGICAL DEVELOPMENT OF THE CONDITIONAL
SAMPLING METHOD. PART I: SENSITIVITY TO STATISTICAL
AND TECHNICAL CHARACTERISTICS

A. K. FOTIADI¹, F. LOHOU¹, A. DRUILHET¹, D. SERÇA¹,*, Y. BRUNET²,
and R. DELMAS¹

¹Laboratoire d'Aerologie, UMR CNRS-UPS: 5560, 14 avenue Edouard Belin, 31400 Toulouse, France; ²INRA - Bioclimatologie, BP 81, 33883 Villenave d'Ornon, France

**Abstract.** The relaxed eddy accumulation (REA), method based on the conditional sampling concept, has received increasing attention over the past few years as it can be used to measure surface fluxes of a wide variety of trace gases for which fast response analysers are not available. In the REA method, a turbulent flux is simply expressed as the product of the standard deviation of vertical wind velocity, the difference between mean scalar concentration in the updrafts and downdrafts and an empirical coefficient, *β* (about 0.63 as based on simulations with a Gaussian distribution, and 0.58 as derived from experimental data). A simulation technique is developed here to evaluate the performance of a ground-based REA system. This analysis uses generated series whose internal structure can be controlled to a large extent. They are stationary and their characteristics are similar to those of physical turbulence. In a first step the influence of some statistical characteristics of vertical velocity and scalar concentration series is investigated. The effect of the third- and fourth-order moments can explain to some degree the difference between calculated and measured *β* values. The impact of a threshold on the vertical velocity is then considered, and its effect on *β* is quantified. The influence of the time lag between *w* and the effective scalar sampling, and the consequences of lowpass filtering of the *w* signal are also investigated. The simulation technique presented in this study can be used to develop elaborate algorithms for near real-time conditional sampling, based on the statistical characteristics of the previous sample.

**Keywords:** Atmosphere-surface exchange, Conditional sampling, Flux measurement, Relaxed eddy accumulation, Turbulent functions.

**1. Introduction**

In the last decade, increasing concern about global climate change and the
enhanced greenhouse effect has required accurate assessment of surface
sources and sinks of chemical compounds (e.g., CO2, CH4, O3, N2O, NOx,
volatile organic compounds, VOC) and determination of the energy balance
over various surfaces. Natural ecosystems can be both sources and sinks of a
great variety of trace gases. Biosphere-atmosphere exchanges of energy and

<footnote>* E-mail:serd@aeroobs-mip.fr</footnote>
mass are controlled by turbulent transfer in the surface boundary layer and the micrometeorological methods developed for measuring energy fluxes can also be applied to trace gas fluxes. The most direct technique to measure turbulent fluxes is the eddy correlation (EC) method by which a vertical flux is determined as the covariance between the fluctuations in the vertical velocity (w) and scalar concentration (X). The EC method requires fast-response sensors (at least 10 Hz) that are not available for many trace gases, thereby limiting its application. However, EC has been successfully used for measuring both energy fluxes (sensible and latent heat) and the fluxes of more conventional trace gases such as CO₂ and CH₄ (Beverland et al., 1996a, b) or VOC (Baker et al., 1999; Pattey et al., 1999; Christensen et al., 2000; Gallagher et al., 2000).

The relaxed eddy accumulation (REA), initially proposed by Businger and Oncley (1990), provides an alternative approach in the absence of fast-response sensors. It involves sampling the air at a constant flow rate and dividing the stream into two different reservoirs, depending on the sign of the vertical wind velocity. At the end of the sampling period the air accumulated in each reservoir is analysed with slow response instruments. The mean scalar flux is then assumed proportional to the concentration difference between the two reservoirs and is given by

$$ \Phi_{wX} = \beta\sigma_w(\overline{X^+} - \overline{X^-}) = \beta\sigma_w\Delta X, \quad (1) $$

where $\sigma_w$ is the standard deviation of the vertical velocity (m s⁻¹) during the sampling period, $\overline{X^+}$ and $\overline{X^-}$ are the mean concentrations of the scalar in updrafts ($w > 0$) and downdrafts ($w < 0$), respectively, and $\beta$ is a dimensionless empirical coefficient. Considering the bi-dimensional distribution of the w and scalar values and making some simple assumptions (normal distribution of w and linearity in the regression of X versus w), Baker et al. (1992) and Baker (2000) estimated the value of $\beta$ as 0.627. Katul et al. (1994, 1996) and Nie et al. (1995) simplified this statistical analysis by taking into account the properties of vertical wind velocity alone and estimated $\beta$ as about 0.63 (0.625 and 0.63, respectively). However, lower values (ranging between 0.56 and 0.60) have been obtained by simulating the REA method from real w and X turbulent time series (Businger and Oncley, 1990; MacPherson and Desjardins, 1991; Baker et al., 1992; Majewski et al., 1993; Pattey et al., 1993; Gao, 1995; Beverland et al., 1996a; Katul et al., 1996). Baker et al. (1992) and Katul et al. (1996) attributed this difference in the $\beta$ values to the departure of the w distribution from a Gaussian one. Katul et al. (1996) and Milne et al. (1999, 2001) discussed the impact of the statistical moments of w and X on $\beta$. For instance Katul et al. (1996), using experimental time series for w and estimating $\beta$ by a model originally proposed by Baker et al. (1992), found that $\beta$ increases with the kurtosis of w. Hence, according to them, the fact that measured $\beta$ are less than 0.62 cannot be attributed to the effect of the w kurtosis. Milne et al. (1999, 2001) investigated the role of the higher-order cross-moments on the $\beta$ values,
using experimental data and a Gram-Charlier two-dimensional probability
distribution analysis. An analytical solution for $\beta$ is presented in terms of two
fourth-order moments: $w'^4$ and $X'w'^3/r$ where $r$ is the correlation coefficient
between the two functions $X$ and $w$. Milne et al. (2001) found a $\beta$ value of 0.557
using experimental values of the fourth-order moments. Other works suggested
that $\beta$ is not significantly influenced by atmospheric stability (Businger and
Oncley, 1990; MacPherson and Desjardins, 1991; Baker et al., 1992; Gao, 1995;
Beverland et al., 1996a,b; Katul et al., 1996) and is relatively insensitive to the
turbulence intensity (Baker et al., 1992). However, Andreas et al. (1998)
showed that the $\beta$ coefficient calculated from momentum fluxes depends on
surface-layer stability expressed in terms of $z/L$, $z$ being the measurement
height, and $L$ the Obukhov length. It varies from 0.56 to 0.58 under unstable
conditions and has a value of 0.63 in neutral conditions. They also found that
the $\beta$ values derived from sensible and latent heat flux depend on stability,
though more weakly than in the previous case. Other specific approaches have
been conducted, showing that $\beta$ could be sensitive to the measurement height
within and above a forest canopy (Gao, 1995).

The conditional sampling technique, as in all other ‘turbulence’ methods,
is limited by different constraints related to homogeneity and stationarity. It
also requires horizontally homogeneous and flat surfaces with a sufficient
fetch to achieve fluxes representative of an homogeneous ecosystem. On top
of this, REA behaves like a ‘black box’, since air sampling, driven by the
vertical velocity signal, is performed in real time. Thus, it is not possible to
further monitor the sampling process in real time. However, if the *w* function
is recorded, its statistical characteristics can be taken into account *a posteriori*
in the flux calculation. On the contrary, the time variation of the scalar
function is not available and only the assumption of similarity with other
scalar fluxes (such as the water vapour flux) can allow one to inspect *a
posteriori* the estimated scalar flux values.

The necessity to quantify the impact of *w* and scalar properties on flux
estimation by REA is therefore obvious. In order to define or test algorithms
of air selection (threshold, signal filtering, etc.), to study the impact of sta-
tistical characteristics of the turbulent functions and to define the application
domain of the method, preliminary studies were made using both generated
times series and experimental data.

Experimental functions allow the influence of atmospheric stability to be
analysed and, through spectral analysis, provide more information than a
relatively simple statistical analysis of turbulent moments. Additionally, using
experimental data, we can go further into the turbulence statistics, either using
functions characterised by a singular coupling of the third- and fourth-order
moments or by studying the chronological succession of turbulent events and
the asymmetry between the length scales of turbulent upward and downward
motions. Simulated turbulence series present two advantages: (i) they can be
perfectly homogeneous and stationary, as compared with experimental
functions and (ii) they allow one to perfectly control and systematically
explore a large range of their statistical characteristics (i.e., skewness and
kurtosis), presumed to have an impact on the REA measurement. Although
simulated series miss the diversity of the physical ones, the way turbulent
series are represented (canonical form or Fourier spectrum) allows their
spectral structure and mainly their energy (variance) to be controlled. They
are therefore well suited to sensitivity analysis of a whole range of factors.

This paper uses simulations with generated series in order to investigate a
wide range of variables. The performed simulations give information about the
impacts of statistical characteristics of the functions or technical choices over the
flux value and then, how the β value should be adapted. As this study intends to
be a practical help to those who make measurements, all the variables tested are
deliberately chosen to be easily calculated. The study has three main objectives:
(i) to investigate the sensitivity of β to the correlation between vertical velocity
and scalar function, (ii) to quantify the impact of statistical characteristics
(third- and fourth-order moments) of w and scalars on a ground-based REA
system and (iii) to test the sensitivity of the method to a threshold on w, lowpass
filtering of the w signal and time lags between w and the scalar concentration.

In a companion paper (Part II), a methodological approach for the analysis
of experimental data will be developed. This analysis is based on a number of
criteria defined so as to allow a quality control of the REA measurements.

## 2. Simulated Turbulent Series

Two kinds of series, coherent and incoherent, were generated to carry out the
simulations. The main difference between the two consists in the presence or
the absence of an internal correlation scale. Incoherent series are easier to
generate and are suitable for studying the sensitivity of the method to various
statistical characteristics of the series. On the contrary, testing of the selection
algorithms (such as filtering or thresholds) requires series presenting an
internal correlation scale. Both kinds of generated series can have a Gaussian
distribution. If non-Gaussian statistics are needed, the series can be further
manipulated so as to meet any pre-defined criterion. As only the results
obtained with coherent series are analysed here, the method to generate
incoherent series is not presented.

### 2.1. SERIES GENERATION

Series are generated from an analytical energy spectrum whose characteris-
tics are very close to those of physical turbulence. More specifically, this
model spectrum has an inertial subrange described by a −5/3 power-law. The
frequency corresponding to the scale of energy production (spectral peak), is pre-defined. It can be shifted to low frequencies in order to simulate the presence of larger-scale processes in the turbulent functions. Application of an inverse Fourier transform allows turbulent series to be generated, with their phases randomly chosen. Coherent series are generated by $(w, X)$ pairs with a given correlation. This is achieved by making their phases dependent on each other in a particular frequency band $\Delta f = f_{max} - f_{min}$, which corresponds to the chosen energetic domain of the cospectra. After the phase function of $w((\varphi_w))$ is randomly chosen, the phase of $X, \varphi_X$, is defined by

$$ \varphi_X = \phi(f) \times \varphi_w + (1 - \phi(f)) \times R, \qquad (2) $$

where $R$ is a frequency random function and $\phi(f)$, a function of the frequency, which is equal to zero, apart from the frequency band $\Delta f$. In the range of $\Delta f$, $\phi(f)$ is symmetric. The correlation degree depends on the width of the common frequency band and its relative position in the frequency domain. The larger the common band is, the higher the correlation between the two series is, provided that the spectral peak belongs to this band. Figure 1a shows two spectrum models with different peak frequencies. Solid and dashed line spectra are used to generate the $w$ and $X$ series respectively, shown in Figure 1b. It is clear here that the dashed line spectrum, whose maximum has been shifted to lower frequencies, provides a series that contains events with more low-frequency energy. The relation between the width of the common band frequency and the correlation coefficient between the two series is illustrated in Figure 1a. The low ($f_{min}$) and high ($f_{max}$) frequency limits of the band are indicated by the vertical dashed line and the bar symbols, respectively. As mentioned above, changing the width of the common band by shifting either the dotted line or the bar symbol allows control of the degree of correlation between the two series. The bar ordinate (y-axis on the right) corresponds to the correlation coefficient of the series. In the case of Figure 1a, the common frequency band has been widened by shifting $f_{max}$ towards the high frequencies, implying an increase of the correlation coefficient from 0.24 to 0.69. The correlation between the two coherent series illustrated in Figure 1b is 0.42.

To sum up, the generated series have similar spectral characteristics to measured functions ones, since the pre-defined spectrum verifies the main characteristics as the observed ones. Concerning the individual statistical moments, their values are very closed to a Gaussian function. At last, the correlation coefficient between the two generated series is under control.

## 2.2. ESTIMATION OF $\beta$

A set of 300 pairs $(w$ and $X)$ of time series, with 4000 points each, was generated. The correlation coefficient varies from 0 to nearly 0.8. For each
<figure><img src="image_2.png" /><figcaption>Figure 1. (a) Example of the energy spectrum used to generate coherent series. The ordinate units are arbitrary. The vertical dashed line indicates the low limit frequency limit of the common band. The abscissa bar symbols show the high frequency limit. The distance between the dashed line and the bar symbol corresponds to the bandwidth (horizontal grey lines). The bar ordinate (y-axis on the right) indicates the correlation between the two generated functions (numbers in brackets). (b) Example of a couple of coherent functions, with a correlation of about 0.42. Both w and X values are normalised.</figcaption></figure>

pair, the mean vertical turbulent flux of a scalar can be expressed as in the EC method. This gives

$$ \Phi_{wX} = \overline{w'X'} = r_{wX}\sigma_w\sigma_X, \qquad (3) $$
<figure><img src="image_3.png" /><figcaption>Figure 2. Concentration difference normalised by the standard deviation versus the correlation between w and scalar concentration. The solid line is the estimated linear regression curve. The inverse of the slope is equal to β.</figcaption></figure>

where w' and X' are the vertical wind velocity and scalar fluctuations respectively, $r_{wX}$ is the correlation coefficient between w and X, and $\sigma_w$ and $\sigma_X$ are the standard deviations of the vertical velocity and scalar concentrations, respectively. By combining (1) and (3), we obtain

$$ \beta = \frac{1}{\Delta X} r_{wX} \sigma_X. \qquad (4) $$

This is the basis for the calculation of individual β values.

The mean value of β is preferably estimated from the slope of the linear regression of $\Delta X/\sigma_X$ against $r_{wX}$, performed over the whole dataset of 300 time series. Figure 2 illustrates the relation between $\Delta X/\sigma_X$ and $r_{wX}$. All data points are equally distributed around the regression line within the whole correlation range. The mean β value obtained through this generation process is $0.627 \pm 0.027$ (0.027 being the mean normalised distance of the points to the regression line). As has already been mentioned in the introduction, identical (Baker et al., 1992; Baker, 2000) or similar (Katul et al., 1994, 1996; Nie et al., 1995) values for β have been reported in other works based on studies limited to Gaussian statistics and first-order cross-moments between variables.

This dataset can then be used to evaluate the dependence of β, as estimated by Equation (3), on a series of parameters such as the statistical moments, a threshold on vertical velocity, the time lag between w and X, etc. In most of the following cases, the variables used to represent the results are normalised by a reference value. It should be noted that, except for our study of the impact of the statistical moments on the REA method (Sections 3.2 and 3.3), the generated series have a Gaussian distribution (the skewness is zero and the kurtosis is 3).
<figure><img src="image_4.png" /><figcaption>*Figure 3.* Values of $\beta$ computed by Equation (3) for each $(w, X)$ couple versus the correlation between $w$ and $X$. The solid line corresponds to the mean value (0.62) and diamond symbols represent the mean $\beta$ value for each range of correlation defined by the vertical lines.</figcaption></figure>

# 3. Sensitivity to Statistical Characteristics of the Turbulent Series

## 3.1. SENSITIVITY OF $\beta$ TO THE $w-X$ CORRELATION

Figure 3 shows the variation of $\beta$ with the $w-X$ correlation. It can be seen that for correlation coefficients larger than about 0.2–0.3 the scatter is fairly low; the mean values of $\beta$ (diamonds), averaged over correlation bands of width 0.2, are nearly constant and close to 0.627. Below 0.2–0.3 the scatter increases significantly.

The accuracy reached in the estimation of a statistical moment depends on how and how rapidly the integral by which it is estimated converges to a limit value. For turbulent series this accuracy increases with the integration time, provided that the series used are stationary. For a given location of the spectral peak, the longer the integration time, the better the parameter converges to its limit value. For perfectly stationary series the integral converges monotonically but usually the convergence is not monotonic: the integral oscillates around the limit value. The same remarks hold for the calculation of $\beta$.

In order to study the convergence of $\beta$ as a function of the integration time, its value is estimated from series of varying length. As the peak in the $w$ spectrum is not always located at the same frequency, the series length ($L$) is normalised by the wavelength ($\lambda$) corresponding to the peak in the $w$ spectrum. This allows comparisons between various generated series to be per-
formed in similar conditions. Indeed, the integration time has to be signifi-
cantly larger than the integral length scale in order to ensure a good enough
statistical representation of the largest eddies.

Figure 4a and b show examples of how $\beta$ converges for a $w-X$ correlation of 0.35 and 0.75, respectively. In each case, 20 series whose normalised length increases up to 15 are plotted. It can be seen that in the case of strong correlation $\beta$ values converge more rapidly and are less scattered around their mean value. Additionally, once the convergence is reached, for $L/\lambda > 10$, the residual noise is less important.

In order to provide more quantitative results, several parameters have
been calculated for six different w-X correlations (0.25, 0.35, 0.45, 0.55, 0.65
and 0.75) (Figure 4c): (i) the standard deviation of β (statistics over 20 series)
for L/λ equal to 2.5 (σ-2.5), 5 (σ-5), 10 (σ-10) and 15 (σ-15); (ii) the
standard deviations of β (σ<sub>conv</sub>) for each of the 20 series, for L/λ varying
between 10 and 15; (iii) the mean value (σ<sub>a-conv</sub> of the 20 σ<sub>conv</sub>) values. In fact,
the first standard deviations quantify the scatter of β around its mean value
and the second ones quantify the residual noise of β once convergence is
reached. Figure 4c shows the variation of these parameters with the w-X
correlation. It is clear that all standard deviations exhibit a logarithmic de-
crease with the correlation. In other words, the error on the flux increases
logarithmically with decreasing correlation. The residual convergence noise
(σ<sub>a-conv</sub>) has a maximum of 0.009 for a correlation coefficient of 0.25. This
value represents an error of only 0.14%, which cannot have an important
effect on the flux. Besides, Figure 4c confirms the fact that the scatter of β
around its mean value increases when L/λ decreases. In the best case for
convergence (L/λ = 15) the scatter in β values (σ - 15) is of order 0.3 for a
w-X correlation of 0.25 (a correlation commonly obtained between w and
scalar physical functions). This means that a non-negligible statistical error of
about 14.5% can be introduced in the flux estimation. When L/λ is lower
than 15, the resulting error increases logarithmically with L/λ.

To sum up, the scatter in $\beta$ increases logarithmically with a decrease in the $w-X$ correlation coefficient and a decrease in the normalised length of the simulated series. Even if the mean $\beta$ value is not affected (Figure 3), the statistical convergence error introduced on a given sample is not negligible.

3.2. SENSITIVITY OF $\beta$ TO THE THIRD-ORDER MOMENT

Physical functions are often skewed due to convective instability in the
atmospheric boundary layer. The logarithmic wind profile in the surface layer
is also a source of asymmetry in turbulent functions. To analyse the possible
influence on $\beta$ of a non-Gaussian probability distribution of $w$, it is desirable
to investigate the relation between $\beta$ and the third- and fourth-order
<figure><img src="image_5.png" /><figcaption>Figure 4. (a) Examples of variation in the $\beta$ convergence as a function of the integration length normalised by the wavelength corresponding to the spectral peak, for a w-X correlation of 0.35. (b) Same as for (a) for a w-X correlation of 0.75. (c) Standard deviations of $\beta$ convergence (in log scale) defined for series of different length, with respect to w-X correlation. Standard deviations $\sigma$-2.5, $\sigma$-5, $\sigma$-10 and $\sigma$-15 correspond to $L/\lambda$ equal to 2.5, 5, 10 and 15, respectively. $\sigma_{conv}$ represents the standard deviation of $\beta$ for each of the 20 series, for $L/\lambda$ in the range 10-15. $\sigma_{a-conv}$ is the mean value of $\sigma_{conv}$.</figcaption></figure>
moments (Figures 5 and 6). These moments are expressed by the dimensionless skewness and kurtosis parameters, respectively. In order to perform this analysis, the probability distribution of the (initially Gaussian) generated series is modified. After the inverse Fourier transform, the series skewness is modified through an exponential operator. A skewed distribution is no longer characterised by a kurtosis of 3. Consequently, when the skewness is modified, the kurtosis is automatically modified. The results presented have been obtained using a large number of generated series.

Kurtosis follows an exponential evolution when the skewness of the simulated series is modified from 0 to 2. The kurtosis does not exceed 5 for a skewness lower than 1.2 and reaches 10 when the skewness is set to 2. The effect of the coupled skewness–kurtosis on $\beta$ is different depending on whether $w$ or the scalar series is considered. Figure 5 illustrates, in the two-dimensional (2D) form, the impact of the skewness (in the range 0 to 3) of both $w$ and the scalar series on $\beta$ (keeping in mind that the kurtosis also changes). If only the skewness and kurtosis of $w$ is taken into account, $\beta$ decreases on average from 0.62 to 0.53. On the contrary, scalar skewness and kurtosis induces a slight increase in $\beta$ from 0.62 to 0.63 on average, for a zero $w$ skewness (and a $w$ kurtosis of 3). In the case where both $w$ and scalar series are asymmetrical, $\beta$ decreases from 0.62 to 0.56. In the 2D diagram the most important variation of $\beta$ is observed along the secondary diagonal [($\text{skew-}\,w=0$, skew-$X=3$), ($\text{skew-}\,w=3$, skew-$X=0$)] where $\beta$ varies from 0.63 to 0.53. Extreme values obtained from the simulations range from 0.46 to 0.72. Considering that the skewness of physical turbulence functions usually varies between 0 and 1, values of the order of 3 can be considered extreme. Indeed, one of the advantages of simulations is that they allow one to explore such extreme values of statistical properties of the functions. For this range of

<figure><img src="image_6.png" /><figcaption>Figure 5. Variation of $\beta$ with the skewness of $w$ (x-axis) and $X$ (y-axis).</figcaption></figure>
values the reduction in $\beta$ is only about 4%, which does not entirely explain the difference between ‘theoretical-simulated’ and ‘experimental’ value of $\beta$.

### 3.3. SENSITIVITY OF $\beta$ TO THE FOURTH-ORDER MOMENT

Contrary to the third-order moments, the fourth-order moment can be modified without changing the skewness value. The impact on $\beta$ of the kurtosis of the vertical velocity and the scalar concentration is presented in a two-dimensional form (Figure 6). The kurtosis varies from 3 to 9, which means that only series with a sharp bell shape are considered. It can be seen that $\beta$ slightly decreases with an increase in the w kurtosis, whereas it increases very little with the scalar kurtosis. As a result, $\beta$ remains practically constant when the coupled effect of w and X kurtosis is considered. The maximum variation of $\beta$ over the whole range of kurtosis [(kurt-w = 9, kurt-X = 3), (kurt-w = 3, kurt-X = 9)] is low and does not exceed 3%. The result concerning the impact of the w kurtosis on $\beta$ is not in agreement with Katul et al. (1996) who found experimentally that $\beta$ increases with the w kurtosis. However, the small variation of $\beta$ illustrated in Figure 6 allows us to conclude, as Katul et al. (1996) did, that the kurtosis of the turbulence functions cannot justify by itself why measured $\beta$ values are lower than 0.62. However, the combined effects of skewness and kurtosis, as shown in Figure 5, imply more important $\beta$ variations. It should be pointed out that the sensitivity of $\beta$ to higher-order cross-moments is not quantified here, although these moments are certainly modified. Up to this point we have studied the sensitivity of the REA method to statistical characteristics of w and scalar functions. In what follows, the effect on the REA method of several choices related either to the selection algorithm (threshold on w signal) or to technical characteristics of a REA system (time lag and w filtering) is investigated.

<figure><img src="image_7.png" /><figcaption>Figure 6. Variation of $\beta$ with the kurtosis of w (x-axis) and X (y-axis).</figcaption></figure>
# 4. Sensitivity to Technical Choices

In the following sections, simulated Gaussian series are used. The aim of the simulations is to quantify the relative variation of $\beta$ according to the technical choices, and consequently to determine to what extent the fluxes could be overestimated or underestimated.

## 4.1. SENSITIVITY OF $\beta$ TO THE VERTICAL VELOCITY THRESHOLD

In order to avoid sampling air for small vertical velocities ($\bar{w} \approx 0$), which is likely to induce background noise, a threshold can be imposed on the $w$ signal, defining a deadband. The latter is generally symmetrical around $w = 0$, and air is not sampled (or is diverted to a third reservoir) within this velocity range. The threshold is usually defined as a fraction of $\sigma_w$. Since scalar concentration is correlated with $w$, the introduction of a threshold results in an increase of the scalar concentration difference between the updraft and downdraft reservoirs. This occurs because only the larger eddies (those with sufficient vertical velocity to exceed the threshold), which move air from further along the vertical concentration gradient and contribute strongly to the scalar flux, are taken into account. Hence, sampled air ‘carries’ more extreme events and leads to larger concentration differences. This allows fluxes to be measured with greater accuracy when the resolution of the analyser is limited.

In Figures 7a and b the impact of the threshold on REA sampling is illustrated for various cases of higher-order moments. Four generated series of $w$ with different statistical characteristics are presented: one with a Gaussian probability distribution, an asymmetrical series with high positive skewness, a series with platykurtic distribution (kurtosis < 3) and a series with leptokurtic distribution (kurtosis > 3). In Figure 7a, the cumulative distribution of each series (or integral curve) is represented from 0% to 100%. Figure 7b shows an enlargement of the central part of Figure 7a ($\sigma_w = \pm 1$, cumulative distribution from 20% to 80%). The vertical lines define the rejection band. In this type of representation, the rejected values correspond to the central domain of the cumulative distribution while the fluctuations considered for the analysis belong to the (positive and negative) tails of the distribution.

Table I summarises for each series the percentage of retained (‘negative’ and ‘positive’ tails) and rejected (central band) fluctuations, for a threshold velocity of $\pm 0.4 \sigma_w$. The proportion of retained points varies with the type of series. For a Gaussian series 68% of the points are conserved and are equally distributed in the two tails. A large (positive) skewness introduces a large asymmetry in the distribution of the points retained: 40% of the values are conserved in the ‘negative’ tail, and only 28% in the ‘positive’ tail. In other words the same percentage as for the Gaussian distribution is retained (68%)
<figure><img src="image_8.png" /><figcaption>Figure 7. (a) Cumulative distribution of four generated series with different statistical characteristics: Gaussian series, a series with large positive skewness, a series with kurtosis > 3 and a series with kurtosis < 3. Cumulative distribution is drawn from 0% to 100%. (b) A zoom on the central part of the integral curves corresponding to the events eliminated by the threshold (represented by the vertical lines).</figcaption></figure>

but updrafts and downdrafts are not equally represented in the statistics. When the kurtosis is lower than 3, the introduction of a threshold eliminates a large number of points (52%) since the latter are grouped around the centre of the distribution. On the contrary, when the kurtosis is larger than 3 the threshold is not as efficient (only 24% of the points are rejected). The proportion of retained data, 76% for kurtosis < 3 and 48% for kurtosis > 3, differs significantly from the Gaussian one (68%).

These examples illustrate the role of a vertical velocity threshold and demonstrate the effect of a symmetrical threshold on the percentage of rejected values, in concert with the effect of higher-order moments of the w series. When the threshold is imposed on a skewed series, the distribution of
<caption>TABLE I</caption>

Percentage of points retained and rejected by a symmetrical threshold of ±0.4σ<sub>w</sub> for each function of Figure 7a.

<table><thead><tr><th rowspan="2"></th><th colspan="3">% of points</th></tr><tr><th>'Negative' tail<br>(w &lt; threshold)</th><th>'Positive' tail<br>(w &gt; threshold)</th><th>Rejected points</th></tr></thead><tbody><tr><td>Gaussian</td><td>34</td><td>34</td><td>32</td></tr><tr><td>Skewness &gt; 0</td><td>40</td><td>28</td><td>32</td></tr><tr><td>Kurtosis &lt; 3</td><td>38</td><td>38</td><td>24</td></tr><tr><td>Kurtosis &gt; 3</td><td>24</td><td>24</td><td>52</td></tr></tbody></table>

the 'positive' and 'negative' events retained is modified, whereas when it is applied to series with a kurtosis different from the Gaussian value, it changes significantly the number of retained points. Considering that experimental w functions are usually skewed or have a kurtosis different from 3, these results suggest that a fixed, symmetrical threshold may not be the most appropriate choice. Therefore, for a more equitable representation of the data points, the threshold should be adapted to higher-order moments (i.e., it should be variable and asymmetrical). Bowling et al. (1999) used a hyperbolic threshold (symmetric and asymmetrical) and found that the concentration difference can increase up to a factor of 2.7, as compared with a conventional REA method. Using an asymmetrical threshold should allow the statistical characteristics of the w function to be taken into account. However, these statistical moments should be known at the time of the sampling itself, which is impossible since the w statistics can only be computed at the end of the integration period. One solution could be to use the skewness and kurtosis of the previous sample, at least in relatively stationary atmospheric conditions.

In order to illustrate the effect of an asymmetrical threshold, we again consider the four previous series. This time, the threshold values are appropriately chosen according to the w statistical moments, in such a way that the retained points are symmetrical in the two distribution tails (Table II), and equal to half the total number of points retained in the Gaussian case (i.e., 34% for a threshold of 0.4σ<sub>w</sub>). Table II shows the values of positive and negative threshold required to achieve this goal for the four series. It can be seen that the thresholds vary from 0.2 to 0.52 σ<sub>w</sub>.

Figure 8 illustrates the effect of the velocity threshold on the normalised concentration difference ΔX/σ<sub>X</sub>. In this case we only used simulated series with Gaussian distributions, and a symmetrical threshold around w = 0. It is obvious from Figure 8, where the threshold is normalised by σ<sub>w</sub>, that the difference in the concentrations measured in the reservoirs increases significantly with the threshold, implying a more precise flux estimation. This
<caption>TABLE II</caption>

Threshold values required to get the same percentage of data points (34%) in each distribution tail.

<table><thead><tr><th></th><th>Negative threshold/σ<sub>w</sub></th><th>Positive threshold/σ<sub>w</sub></th></tr></thead><tbody><tr><td>Gaussian</td><td>0.40</td><td>0.40</td></tr><tr><td>Skewness > 0</td><td>0.50</td><td>0.21</td></tr><tr><td>Kurtosis < 3</td><td>0.52</td><td>0.52</td></tr><tr><td>Kurtosis > 3</td><td>0.20</td><td>0.20</td></tr></tbody></table>

<figure><img src="image_9.png" /></figure>

Figure 8. Impact of the threshold on the mean concentration of updrafts (diamonds) and downdrafts (circles). Concentration and threshold are normalised by the standard deviation of the scalar and w series, respectively. Note that since X is centered, downdrafts are associated with negative X.

increase must be compensated by a decrease in $\beta$, illustrated in Figure 9 for a range of thresholds from 0 to $1 \sigma_w$, where $\beta$ is normalised by its reference value (i.e. with no threshold). The resulting curve can be described by the following exponential function:

$$ \frac{\beta}{\beta_0} = e^{-0.67(w_0/\sigma_w)} \qquad (6) $$

with $w_0$ = threshold. Similar results were reported by Businger and Oncley (1990), Pattey et al. (1993) and Katul et al. (1996) for experimental temper-
<figure><img src="image_10.png" /><figcaption>Figure 9. Sensitivity of $β$ to the threshold. Diamond symbols represent our results derived from simulations with generated Gaussian series. Square and circle symbols show the $β$ values resulting from the models of Businger and Oncley (1990) and Pattey et al. (1993), respectively. $β$ is normalised by its reference value, corresponding to zero threshold.</figcaption></figure>

ature, water vapour and CO$_2$ data. Some of their results are also illustrated in Figure 9. For threshold values lower than 0.3 $σ_w$ all models give similar estimates for $β$. For instance, for threshold of 0.1 $σ_w$ and 0.2 $σ_w$ typically used in operational REA systems, $β$ has to be reduced by about 8% and 15%, respectively. Deviations between the three models for the high values of the threshold can be explained by the statistical characteristics of the *w* series (different departures from the Gaussian distribution).

From a technical point of view, it has been suggested that a threshold can be used to extend the life of the valves by limiting their switching frequency (Pattey et al., 1993). Indeed, the use of a threshold leads to a decrease in the percentage of the total time of updraft and downdraft sampling (Figure 10). In order to investigate the role of the threshold on valve operation, the number of valve switchings is counted for different thresholds. Six types of switching are defined. Firstly, two switchings called ‘direct’ represent transitions from positive to negative *w* values and vice versa, and secondly four ‘indirect’ switchings correspond to transitions to and from the deadband. As the simulated series have, by construction, perfectly Gaussian distributions, indirect transitions are equally probable and their curves are identical. In Figure 11 all indirect transitions are therefore represented by just one curve
(square symbols). The same holds for direct transitions (circles). The line with
diamonds represents the total number of transitions. As expected, Figure 11
shows that the number of direct transitions decreases monotonically when
the threshold increases, whereas the number of indirect transitions and the
total number of transitions first increase up to a threshold of 0.5 σ<sub>w</sub> and 0.4
σ<sub>w</sub>, respectively, then decrease. However, even for a large threshold (1.0 σ<sub>w</sub>),
the total number of transitions remains greater than that corresponding to a
zero threshold. A limitation in the number of transitions could be achieved
for higher values of the threshold (outside the range represented in Figure 11)
but the number of rejected values would become too important. Therefore,
the use of a threshold does not really contribute to optimise, in terms of
decrease, the switching frequency of the valves.

To conclude, it should be pointed out that using a threshold allows flux
measurements to be performed when the concentration difference is small rel-
ative to the accuracy of the analyser. However, it does not reduce the total
number of valve switchings. Also, there is potential advantage in using an
asymmetrical threshold, because when a symmetrical threshold is applied to
skewed and flattened series the risk of reducing significantly the effective sam-
pling duration (by eliminating a large number of air samples) is not negligible.

4.2. SENSITIVITY OF β TO A TIME LAG BETWEEN w AND X

Another source of uncertainty that should be considered with the conditional
sampling method is the delay in the time response of the REA system to
changes in the sign of w and its consequences on the effective conditional

<figure><img src="image_11.png" /><figcaption>*Figure 10.* Impact of the threshold on the total percentage of time for positive sampling (squares), negative sampling (triangles) and no sampling (circles).</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>Figure 11. Frequency of valve switchings versus the normalised threshold. Circle and square symbols represent the direct and indirect transitions, respectively. The diamonds describe the evolution of the total number of transitions.</figcaption></figure>

sampling of air. There are three main causes for this delay: (i) the time lag between the change in vertical velocity sign and the detection of this change by the processing system, (ii) the time lag between the detection of this change and the subsequent actuation of the sampling valves and (iii) the length of the inlet tubing, which may have an important contribution to the time lag effect as it takes a finite amount of time for the air sample to travel through the sample line until it reaches the valves. The total time delay induces an underestimation of the concentration difference and, consequently, of the estimated flux since a fraction of the updraft and downdraft air can be routed to the wrong reservoir. This delay can vary from a few milliseconds to a few seconds and obviously depends on the characteristics of each system. The underestimation of the concentration difference due to the time lag is also a function of the valve switching frequency (the mean frequency with which *w* changes sign) or eddy-reversal frequency as defined by Baker et al. (1992). In its turn, mean switching frequency directly depends on the turbulence conditions.

Time lags related to the time response of the instruments also exist in eddy correlation systems. In this case the time lag is easily estimated through the calculation of the cross-correlation between *w* and *X*, as the time lag for which the cross-correlation shows a maximum. The flux is then calculated from the lagged functions. This is not possible for REA measurements since the sampled air is stored and analysed later. However, the underestimation of
$\Delta X$ can be balanced by adjusting the $\beta$ coefficient. In order to determine the appropriate adjustment of $\beta$, simulations were made for a time lag ranging from 0 to 6.4 s, so as to include the internal scale of correlation.

The dependence of some REA parameters on the time lag is illustrated in Figure 12. The square curve illustrates the variation of the mean ‘positive’ concentration (scalar concentration in the ‘positive’ reservoir corresponding to updrafts) with the delay. Note that for series with a Gaussian distribution the mean ‘negative’ concentration curve is exactly symmetrical. Consequently, the positive curve represents the concentration difference between the two reservoirs, divided by two. Circle and diamond curves represent the variations in $\beta$ and the $w-X$ lagged correlation, respectively. These results represent averages over 20 simulated series, and the error bars on the $\beta$ curve represent the corresponding standard deviations. The correlation between $w$ and $X$ pairs is of the order of 0.5. In order to represent measurements carried out at a 6-m height, a 30-m spectral length scale is assumed. With a typical mean horizontal wind of $5 \text{ m s}^{-1}$ the peak frequency is then 0.13 Hz. Figure 12 shows that the REA parameters are strongly influenced by the time lag. Indeed, as the time lag increases the concentration difference decreases, along with a concomitant decrease in the $w-X$ correlation. Such variations result in a large and rapid increase in $\beta$ (by 3%, 11% and 26% for a time lag of 0.2, 0.4 and 0.6 s respectively). $\beta$ values become very large for lags greater than 1 s. For a long enough delay (4.4 s) the selection sign is reversed, yielding to a flux

<figure><img src="image_13.png" /><figcaption>Figure 12. Sensitivity of REA parameters to the time lag (in s): correlation between $w$ and $X$ (diamonds), $\beta$ (circles) and mean concentration of updrafts (square). The error bars represent the standard deviation. Sampling rate, 10 Hz.</figcaption></figure>
of opposite direction. The correlation also changes sign and a point of dis-
continuity appears on $\beta$ (not shown here), after which it varies almost ran-
domly. The scatter in the $\beta$ values also increases with the time lag and
becomes infinite when the correlation changes sign. It has to be pointed out
that the percentages in the $\beta$ increase given above cannot be considered as a
general result because they strongly depend on the correlation degree be-
tween $w$ and scalar series. However, this example allows us to illustrate the
order of magnitude of the uncertainty introduced by the time lag in flux
calculation.

4.3. SENSITIVITY OF $\beta$ TO LOWPASS FILTERING

Real-time lowpass filtering is a common practice that is applied to smooth
out fluctuations of the vertical velocity signal (Baldocchi et al., 1988; Beier,
1991; Delon et al., 2000) and thus to increase the signal-to-noise ratio. This
ratio may be low due to the presence of noise, resulting mainly from technical
and physical causes: (i) the instruments and acquisition systems used to
measure turbulent fluxes may add digital noise (aliasing) to the $w$ signal; (ii)
in some cases, like nocturnal stable conditions, high-frequency fluctuations
are relatively more important in the turbulence spectra. Lowpass filtering
also allows one to reduce the switching frequency of the valves. Filtering $w$ in
real-time implies using an asymmetrical filter (in time) as only the past history
of the function is known, and introduces a phase shift in the $w$ function which
in turn alters the correlation between $w$ and the scalar concentration, thereby
affecting the calculated flux.

In order to quantify the effect of low-pass filtering on flux calculation,
simulations were carried out using a binomial recursive filter with a variable
cut-off frequency ranging from 5 to 0.4 Hz. In Figure 13, three variables are
represented against the logarithm of the cut-off frequency: the correlation
between $w$ and $X$ (line with diamonds), the $\beta$ coefficient (line with circles) and
the signal delay introduced by filtering $w$, expressed in seconds (line with
squares). Here again the results represent averages over 20 series, with the
corresponding standard deviations shown as error bars on the $\beta$ curve. It can
be seen that the effect of lowpass filtering is similar to that of an increasing
time lag. It introduces a delay in $w$ signal detection and an attenuation of the
signal amplitude, which increases with a reduction in the cut-off frequency.
Due to the delay of the $w$ function, the correlation between $w$ and $X$ de-
creases, and $\beta$ values and their standard deviation rapidly increase as the cut-
off frequency decreases. It increases by 1% and 7% for a cut-off frequency of
5 and 2.5 Hz, respectively. This increase becomes huge for a cut-off frequency
lower that 1 Hz.
<figure><img src="image_14.png" /><figcaption>Figure 13. Sensitivity of $\beta$ (circles), w-X correlation (diamonds) and w signal delay (squares) to lowpass filtering. The x-axis represents the logarithm of the cut-off frequency.</figcaption></figure>

**5. Conclusions**

This study focuses first on the sensitivity of the REA method to statistical characteristics of the turbulent vertical velocity and scalar series, and second on a number of technical choices associated with the conditional sampling algorithm. This analysis of the REA methodology was performed using artificial, generated series.

A first general conclusion deduced from this study concerns the usefulness of simulations made with generated series. The technique described here appears as a powerful tool for methodological development and validation of new techniques and algorithms aimed at measuring trace gas fluxes. The spectral characteristics of generated series are similar to those of experimental functions and easy to deal with. Their statistical characteristics (skewness and kurtosis) are originally close to the Gaussian ones but can be widely modified. Apart from the obvious advantage of this kind of simulation (the analysis and understanding of the functioning mode of the REA method), such series have the advantage, compared to measured data, to allow us to investigate the behaviour of the method to extreme values of the parameters. As such, they provide a nice complement to experimental functions.

Our simulations with Gaussian series provide a mean value of 0.627 for the $\beta$ coefficient, with no dependence on the correlation between the vertical velocity and the scalar concentration. The large scatter obtained for the low correlation values can be attributed to the relatively slow convergence of the
integral leading to the estimate of the mean $\beta$ value. The convergence towards the mean value is a function of the w-X correlation and also depends on the integration time (length of the time series, expressed as a function of $L/\lambda$). The convergence is usually reached for $L/\lambda > 10$. After this limit, the residual noise depends on the w-X correlation value.

We summarise in Table III the results quantifying the uncertainties introduced in the flux calculation (and the resulting required adjustment of $\beta$) by (i) the statistical characteristics (third- and fourth-order moments) of w and scalar functions, (ii) the vertical velocity threshold, (iii) the time lag between w and the effective sampling and (iv) the lowpass filter applied on the w signal. The values in italics relate either to cases that can be encountered in experimental data (skewness and kurtosis), or to cases used in practice (threshold and filtering). The main conclusions are as follows.

* Positive third- and fourth-order moments of w induce a decrease in $\beta$ by about 3–4% on average, for a skewness around 1. A positive scalar skewness and kurtosis cause a slight increase. Their combined effect results in a decrease of $\beta$.

* The fourth-order moment alone has no major effect on $\beta$. A slight decrease with an increase in w kurtosis and a slight increase with an increase in the scalar kurtosis have been found. Kurtosis cannot explain why the experimentally derived values of $\beta$ are consistently smaller than 0.62.

* The introduction of a threshold on the vertical velocity enhances the method sensitivity by increasing the concentration difference between updrafts and downdrafts. An exponential relationship describes the impact of the threshold on $\beta$, in agreement with previous studies. The impact on $\beta$ of the most common thresholds used in the field has been quantified (Table II). On the contrary, the switching frequency of the valves first increases with the threshold, then decreases. The impact of a symmetrical threshold on the number of rejected air samples, when the turbulent series have non-Gaussian skewness or kurtosis, is also demonstrated. For a symmetric threshold of $\pm0.4\sigma_w$, the percentage of rejected points is 24 for kurtosis lower than 3, 52 for kurtosis higher than 3, and 32 for a Gaussian or skewed series. Another test shows how asymmetric the threshold must be to keep the same percentage of data points (34%) in each distribution tail. For a positive skewness, the normalised negative and positive thresholds are equal to 0.5 and 0.2 respectively.

* The time lag (between changes in w sign and effective air sampling) and lowpass filtering have similar effects. The correlation between w and X decreases and $\beta$ values increase very rapidly. For instance for a delay of 0.2 s, a typical value for a ground-based REA system, $\beta$ increases by 3% relative to its value at zero delay (0.62).
<caption>TABLE III</caption>

Errors induced in the flux estimate by the skewness and kurtosis of *w* and *X* functions, the threshold, time lag and lowpass filtering.

<table><thead><tr><th></th><th>Flux error (%)</th></tr></thead><tbody><tr><td>Skewness of <i>w</i></td><td></td></tr><tr><td>1</td><td>-3.2</td></tr><tr><td>3</td><td>-14.5</td></tr><tr><td>Skewness of <i>X</i></td><td></td></tr><tr><td>1</td><td>+0.5</td></tr><tr><td>3</td><td>+1.8</td></tr><tr><td>Skewness of <i>w</i> and <i>X</i></td><td></td></tr><tr><td>1</td><td>-2.1</td></tr><tr><td>3</td><td>-9.8</td></tr><tr><td>Kurtosis of <i>w</i></td><td></td></tr><tr><td>3</td><td>0</td></tr><tr><td>9</td><td>-2.3</td></tr><tr><td>Kurtosis of <i>X</i></td><td></td></tr><tr><td>3</td><td>0</td></tr><tr><td>9</td><td>+2.6</td></tr><tr><td>Kurtosis of <i>w</i> and <i>X</i></td><td></td></tr><tr><td>3</td><td>0</td></tr><tr><td>9</td><td>+0.3</td></tr><tr><td>Threshold (fraction of σ<sub>w</sub>)</td><td></td></tr><tr><td>0.1</td><td>+8</td></tr><tr><td>0.2</td><td>+15</td></tr><tr><td>0.5</td><td>+30</td></tr><tr><td>Time lag (s)</td><td></td></tr><tr><td>0.2</td><td>-3</td></tr><tr><td>0.4</td><td>-11</td></tr><tr><td>0.6</td><td>-26</td></tr><tr><td>Cut-off filtering frequency (Hz)</td><td></td></tr><tr><td>5.0</td><td>-1</td></tr><tr><td>2.5</td><td>-7</td></tr></tbody></table>

In practice, application of these conclusions to a true REA system, for instance by using a non-symmetric threshold, would require real-time calculation of the *w* statistical characteristics, which is not feasible. The solution suggested here is to use the statistics inferred from the previous sampling sequence, which should not be detrimental under usual atmospheric conditions. It is always possible, and this should be done systematically, to check this assumption afterwards and, whenever necessary, discard a par-
ticular run if the statistical characteristics of the turbulent time series turn out
to differ by a large extent from those of the previous run.

As mentioned above, generated series allow the separate influence of
various factors to be analysed, and a large range of conditions to be ac-
counted for. Physical functions cannot be used for the same purpose. They
are more complex than artificial ones; in particular, simultaneous time series
such as vertical velocity and scalar series exhibit complex interdependency
patterns that cannot be reproduced with generated series. In other words,
generated series do not replace experimental data but should be considered as
a powerful complement to them. In a companion paper (Fotiadi et al., 2005)
we will go back to the physical world and use a large set of experimental
turbulent time series to further analyse the REA methodology.

References

Andreas, E. L., Reginald, J. H., James, R. G., Douglas, I. M., William, D. O., and Sarma, A. D.: 1998, ‘Stability Dependence of the Eddy-Accumulation Coefficients for Momentum and Scalars’, *Boundary-Layer Meteorol.* **86**, 409–420.

Baker, B., Guenther, A., Greenberg, J., Goldstein, A., and Fall, R.: 1999, ‘Canopy Fluxes of 2-Methyl-3-buten-2-ol over a Ponderosa Pine Forest by Relaxed Eddy Accumulation: Field Data and Model Comparison’, *J. Geophys. Res.* **104**, 26.107–26.114.

Baker, J.: 2000, ‘Conditional Sampling Revisited’, *Agric. For. Meteorol.* **104**, 59–65.

Baker, J. M., Norman, J. M., and Bland, W. L.: 1992, ‘Field-Scale Application of Flux Measurement by Conditional Sampling’, *Agric. For. Meteorol.* **62**, 31–52.

Baldocchi, D. D., Hicks, B. B., and Mayers, T. P.: 1988, ‘Measuring Biosphere–Atmosphere Exchanges of Biologically Related Gases with Micrometeorological Method’, *Ecology* **69**, 1331–1340.

Beier, N.: 1991, ‘Measuring Fluxes of Chemical Components by Eddy Accumulation’, in 7th *Symposium on Meteorological Observations and Instrumentation*, American Meteorological Society, Boston, MA, pp. 1–5.

Beverland, I. J., Milne, R., Boissard, R., Ónéill, C., Moncrieff, J. B., and Hewitt, C. N.: 1996a, ‘Measurement of Carbon Dioxide and Hydrocarbon Fluxes from a Sitka Spruce Forest Using Micrometeorological Techniques’, *J. Geophys. Res.* **101**, 22.807–22.815.

Beverland, I. J., Moncrieff, J. B., Ónéill, C., Hargreaves, K. J., and Milne, R.: 1996b, ‘Measurement of Methane and Carbon Dioxide Fluxes from Peatland Ecosystems by the Conditional-Sampling Technique’, *Quart. J. Roy. Meteorol. Soc.* **122**, 819–838.

Bowling, D. R., Delany, A. C., Turnipseed, A. A., Baldocchi, D. D., and Monson, R. K.: 1999, ‘Modification of the Relaxed Eddy Accumulation Technique to Maximize Measured Scalar Ratio Differences in Updrafts and Downdrafts’, *J. Geophys. Res.* **104**, 9.121–9.133.

Businger, J. A. and Oncley, S. P.: 1990, ‘Flux Measurement with Conditional Sampling’, *J. Atmos. Oceanic Tech.* **7**, 349–352.

Christensen, C. S., Hummelshøj, P., Jensen, N. O., Larsen, B., Lohse, C., Pilegaard, K., and Skov, H.: 2000, ‘Determination of the Terpene Flux from Orange Species and Norway Spruce by Relaxed Eddy Accumulation’, *Atmos. Environ.* **34**, 3057–3067.

Delon, C., Druilhet, A., Delmas, R., and Greenberg, J.: 2000, ‘Aircraft Assessment of Trace Compound Fluxes in the Atmosphere with Relaxed Eddy Accumulation: Sensitivity to the Conditions of Selection’, *J. Geophys. Res.* **105**, 20.461–20.472.
Fotiadi, A. K., Lohou, F., Druilet, A., Serça, D., Said, F., Laville, P., and Brut, A.: 2005, 'Methodological Development of the Conditional Sampling Method. Part II: Quality Control Criteria of Relaxed Eddy Accumulation Flux Measurement'. *Boundary-Layer Meteorol*, in press.

Gallagher, M. W., Clayborough, R., Beswick, K. M., Hewitt, C. N., Owen, S., Moncrieff, J., and Pilegaard, K.: 2000, 'Assessment of a Relaxed Eddy Accumulation for Measurements of Fluxes of Biogenic Volatile Organic Compounds: Study over Arable Crops and a Mature Beech Forest', *Atmos. Environ.* **34**, 2887–2899.

Gao, W.: 1995, 'The Vertical Change of Coefficient *b* Used in the Relaxed Eddy Accumulation Method for Flux Measurement above and within a Forest Canopy', *Atmos. Environ.* **29**, 2339–2347.

Katul, G., Albertson, J., Parlange, M., Chu, C. R., and Striker, H.: 1994, 'Conditional Sampling, Bursting, and the Intermittent Structure of Sensible Heat Flux', *J. Geophys. Res.* **99**, 22.869–22.876.

Katul, G., Finkelstein, P. L., Clarke, J. F., and Ellestad, T. G.: 1996, 'An Investigation of Conditional Sampling Method Used to Estimate Fluxes of Active, Reactive, and Passive Scalars', *J. Appl. Meteorol. 35*, 1835–1845.

MacPherson, J. I. and Desjardins, R. L.: 1991, 'Airborne Tests of Flux Measurement by the Relaxed Eddy Accumulation Technique', in *Proceeding of the Seventh AMS Symposium on Meteorological Observations and Instrumentation*, New Orleans, LA, January 14–18, 1991, pp. 6–11.

Majewski, M., Desjardins, R. L., Rochette, P., Pattey, E., Selber, J., and Glotfelty, D.: 1993, 'Field Comparison of an Eddy Accumulation and an Aerodynamic-Gradient System for Measuring Pesticide Volatilization Fluxes', *Environ. Sci. Technol.* **27**, 121–128.

Milne, R., Beverland, I. J., Hargreaves, K., and Moncrieff, J. B.: 1999, 'Variation of the $\beta$ Coefficient in the Relaxed Eddy Accumulation Method', *Boundary-Layer Meteorol.* **93**, 211–225.

Milne, R., Mennim, A., and Hargreaves, K.: 2001, 'The Value of the $\beta$ Coefficient in the Relaxed Eddy Accumulation Method in Terms of Fourth-Order Moments', *Boundary-Layer Meteorol.* **101**, 359–373.

Nie, D., Kleindienst, T. E., Arnts, R. R., and Sickles, J. E.: 1995, 'The Design and Testing of a Relaxed Eddy Accumulation System', *J. Geophys. Res.* **100**, 11.415–11.423.

Pattey, E., Desjardins, R. L., and Rochette, P.: 1993, 'Accuracy of the Relaxed Eddy Accumulation Technique, Evaluated Using CO₂ Flux Measurements', *Boundary-Layer Meteorol.* **66**, 341–355.

Pattey, E., Desjardins, R. L., Westberg, H., Lamb, B., and Zhu, T.: 1999, 'Measurement of Isoprene Emissions over a Black Spruce Stand Using a Tower-Based Relaxed Eddy-Accumulation System', *J. Appl. Meteorol.* **38**, 870–877.

Rinne, H. J. I., Guenther, A. B., Warneke, G., De Gouw, J. A., and Luxembourg, S. L.: 2001, 'Disjunct Eddy Covariance Technique for Trace Gas Flux Measurements', *Geophys. Res Lett.* **28**, 3139–3142.