- Email: [email protected]
The impact of sampling rate on eddy-covariance flux estimates
Agricultural and Forest Meteorology 109 (2001) 39–45
The impact of sampling rate on eddy-covariance flux estimates F.C. Bosveld a,∗ , A.C.M. Beljaars b a
Royal Netherlands Meteorological Institute, P.O. Box 201 3730 De Bilt, AE, The Netherlands b European Centre for Medium Range Weather Forecast, Reading, UK
Received 18 September 2000; received in revised form 9 May 2001; accepted 15 May 2001
Abstract It is shown that sampling rate has no influence on the expected value of an eddy-covariance flux estimate. The uncertainty in the flux estimate, however, increases with decreasing sampling rate. This sampling induced uncertainty is compared with the natural uncertainty in the flux estimate due to the stochastic nature of turbulence. It is shown that the sampling induced uncertainty becomes important when the sampling interval is larger than the integral time scale of the flux time series. Theoretical expressions are evaluated by using 9 days of turbulence data taken over a coniferous forest. For sampling intervals smaller than the integral time scale the sampling induced uncertainties are underestimated by the model because the actual auto-correlation function for the flux time series deviates from the exponential form assumed in the model. © 2001 Elsevier Science B.V. All rights reserved. Keywords: Scalar flux; Eddy-covariance; Integral time scale; Turbulence statistics
1. Introduction Vertical fluxes between the earth surface and the atmosphere can be measured directly in the atmospheric surface layer by means of the eddy-covariance method, provided that fast response sensors are available for the determination of the vertical wind speed and the transported quantity under consideration. Over the years significant improvements have been made in correction procedures because of limitations of the instruments and on the minimal requirements of instruments to give reliable flux estimates (Moore, 1986; ∗ Corresponding author. Tel.: +31-30-2206-911; fax: +31-30-2210-407. E-mail address: [email protected] (F.C. Bosveld).
Moncrieff et al., 1997; Massman, 2000). In these studies no attention is given to the subject of sampling rate requirements. The sampling rate required in an experiment depends on the purpose of the study. If one is interested in spectral analysis of the high frequency turbulent regime one is likely to choose the sampling rate related to the time constant of the measuring system. For example, when using a sonic anemometer with a path length of d at a wind speed U one is likely to use a sample rate faster than U/d. Another example is the shift of samples in time to account for delays in measuring tubes (i.e. trace gas analysers). One then would like to sample at a rate with which an accurate shift can be obtained. In such cases required sample rates are typically between 10 and 50 Hz. Note that for spectral purposes low-pass filtering at the Nyquist
0168-1923/01/$ – see front matter © 2001 Elsevier Science B.V. All rights reserved. PII: S 0 1 6 8 - 1 9 2 3 ( 0 1 ) 0 0 2 5 7 - X
40
F.C. Bosveld, A.C.M. Beljaars / Agricultural and Forest Meteorology 109 (2001) 39–45
frequency is desirable to avoid aliasing, whereas for flux estimates this filtering is undesirable due to loss of high frequency flux contribution. For monitoring of vertical fluxes one is often only interested in time averaged fluxes and much lower sampling rates may be chosen. Bosveld and Bouten (2001), for example, use a sampling rate of 1 Hz for eddy-covariance measurements taken over an 18 m high coniferous forest at a height of 30 m. Here we investigate the impact of sampling rate on time averaged flux estimates. Both the impact on the average flux and on the uncertainty of the flux estimate will be quantified. It seems that Haugen (1978) was the first to treat the subject of sampling requirements. Here we follow the treatment of Lenschow et al. (1994) who give a fairly complete analysis. We modify their analysis to avoid interference with the impact of the finite length of the time series, which they treat simultaneously. The latter phenomenon is related to uncertainties in the low-frequency contribution to the flux and is of no concern here. By treating sampling rate separately the nature of its impact on eddy-covariance fluxes will become more clear. The theoretical results are evaluated with turbulent measurements taken at the same location and height as those of Bosveld and Bouten (2001), but with a sampling rate of 10 Hz. Time series with lower sampling rates are constructed from the 10 Hz time series by retaining samples at regular intervals. Fluxes from this constructed time series are then compared with fluxes derived from the original sampled time series. In the following we will assume that spectral losses due to finite response characteristics of the measuring system can be ignored.
First, we derive the expected flux and show that an unbiased estimate is obtained independent of the sampling rate. For the cross product we have {wj (t)sj (t)} = F(t) = F
where F is the flux. The constantness of F follows from the stationarity of the time series over the period T. Each simultaneous measurement of w and s results in an unbiased estimate of the flux and thus the flux calculated from a sampled time series will result also in an unbiased estimate independent of the sampling rate. Note that in practical applications the average of wj (t) and sj (t) over the period T are subtracted from the time series. This results in an extra term in the expression for the flux, representing the loss in the co-variance due to the omission of time scales T and longer. It can be shown that this low frequency term is almost independent of sampling rate, provided that the sampling interval is much smaller than T. In this study, we are only interested in a comparison between time series with different sampling rates and we can omit this term in the analysis. Now, we turn to the relation between the sampling rate and the variance of the flux estimate. For a realisation j an estimate of the flux over the period T is obtained by 1 T FTj = wj (t)sj (t) dt (3) T 0 The expected value of FTj is of course the flux F. Let us define the flux deviation time series as fj (t) = wj (t)sj (t) − F
2. Theory Let wj (t) and sj (t) be time series for t = 0 to T of the vertical wind speed and a scalar quantity, respectively. The time series are supposed to be stationary stochastic processes drawn from an ensemble. The index j indicates the jth realisation of the ensemble. Furthermore, we assume that the expected values (defined as the ensemble average and represented by {·}) of w and s are 0, thus w and s represent the fluctuations. In formula {wj (t)} = 0,
{sj (t)} = 0
(1)
(2)
(4)
This is a stochastic function with an expected value 0 at each instant of time. We now construct a sampled time series with a sampling interval ∆. The averaging period T is divided in N = T /∆ intervals. Samples are taken halfway the sample intervals. This is not a principle point, but this form gives a more elegant transition to the integral formulation. The deviation from F over a period T for the jth realisation of the ensemble when sampled with a time interval ∆ is fjT (∆)
N−1 1 1 = fj n∆ + ∆ N 2 n=0
(5)
F.C. Bosveld, A.C.M. Beljaars / Agricultural and Forest Meteorology 109 (2001) 39–45
41
In the limit for ∆ approaching 0, this transforms into the integral over the period T of the continuous function. Here we are interested in the ensemble variance VfT (∆) of fjT (∆). To this end we write VfT (∆) = {(fjT (∆))2 } N−1 1 1 1 = 2 fj n∆+ ∆ fj m∆+ ∆ 2 2 N n,m=0
(6) This can be written in terms of the auto-correlation function ρff and the variance µf of the stochastic function f N−1 µf ρff ([n − m]∆) VfT (∆) = 2 N
(7)
n,m=0
and with k = n − m and T = N ∆, we obtain VfT (∆) =
µf T
N
ρff (k∆) 1 −
k=−N
|k∆| ∆ T
(8)
and thus, if T is much larger than the time scale on which ρff approaches 0 we can write VfT (∆) =
∞ µf ρff (k∆)∆ T
(9)
k=−∞
By letting ∆ approaching 0 we obtain Eq. (46) of Lenschow et al. (1994). µf ∞ 2µf Jf T Vf (0) = ρff (t) dt = (10) T −∞ T where Jf is the integral time scale by definition. We call this the natural variance of the flux, i.e. the variance induced by the stochastic nature of the turbulence. The sampling induced variance WfT (∆) is defined as the difference between the total variance (Eq. (9)) and the natural variance (Eq. (10)): WfT (∆) = VfT (∆) − VfT (0)
(11)
The summation of Eq. (9) is represented in Fig. 1 for an exponential auto-correlation function with integral time scale Jf and for three different sampling intervals, i.e. ∆ = Jf , Jf /2 and 0. The figure shows that the increase in variance of the sampled case
Fig. 1. The continuous (∆ = 0) and finite sample interval (∆ = Jf and Jf /2) auto-correlation function, ρff , as function of delay time scaled on the integral time scale (t/Jf ), representing the summation of Eq. (9).
compared to the continuous case is mainly caused by the k = 0 term. The other terms in the sampled case are good approximations for the integral of the continuous case in the respective intervals. In fact, it can be shown that for an exponential auto-correlation function the k = 0 terms contribute negatively to WfT (∆) for at most 13%. This maximum contribution is reached when ∆ is of the order of Jf . The magnitude of the k = 0 term decreases with decreasing sampling interval. The figure also shows that, it is not the integral time scale that determines sampling induced variance, but the behaviour of the correlation function at delay times between 0 and ∆/2. For an exponential auto-correlation function with integral time scale Jf an exact form is retrieved by Lenschow et al. (1994), Eq. (58) 2µf Jf ∆ ∆ VfT (∆) = (12) coth T 2Jf 2Jf From this equation we derive that for an exponential auto-correlation function and ∆ equal to Jf the sampling induced variance WfT (∆) is 10% of the natural variance.
3. Experimental evaluation To test the theory on the influence of sampling rate on the mean flux and the variance in the flux,
42
F.C. Bosveld, A.C.M. Beljaars / Agricultural and Forest Meteorology 109 (2001) 39–45
we used a series of 9 days of 10 Hz turbulence data of the vertical wind speed and the sonic temperature, both measured with a Kaijo Denki DAT-300 sonic anemometer–thermometer system. The data were obtained during August 1996 at the Speulderbos research site, a coniferous forest site near the village of Garderen, The Netherlands. The geographical co-ordinates are 52◦ 15 N, 5◦ 41 E and 52 m altitude. Site and instrumentation are described in Bosveld and Bouten (2001). Eddy-covariance measurements were performed at 31 m above the forest floor. With an average tree height at that time of 22 m and a thinned stand the displacement height is estimated at 14 m, resulting in an effective measuring height of 17 m. This value is 2 m smaller than the effective measuring height during the 1989 experiment described in Bosveld and Bouten (2001). Six time series with different sample intervals were obtained by omitting observations. The sample intervals created were 0.1, 0.2, 0.4, 1.0, 2.0 and 4.0 s. The data were processed over 10 min. Sensible heat flux was derived by calculating the covariance’s between the vertical wind speed and the temperature using the six time series with different sample intervals. Fig. 2 shows a time series of heat fluxes for one day for sampling intervals 0.1, 1.0, and 4.0 s. A clear increase in scatter is observed for the larger sampling intervals. To arrive at modelled flux variances from Eq. (12) the variance µf of the fluctuating flux time series and the integral time scale Jf were estimated for each 10 min interval of the 9 daytime series.
Fig. 3. Auto-correlation function of the flux fluctuating time series ρff and of the vertical wind speed ρ ww based on a time period with wind speed of 5 m s−1 and a sensible heat flux of 280 W m−2 . Also shown is an exponential decay function with time constant 1.5 s.
Wyngaard (1973), based on the Kanses experiment, related µf to the flux F. He found µf ≈ 4 F2 for daytime cases and µf ≈ 10 F2 for nighttime cases. These relations were tested for the current data set and confirmed. Together with Eq. (10) we derive for a wind speed of 3.5 m s−1 during daytime a typical uncertainty of 20% in a 10 min averaged flux. Jf was estimated by calculating the auto-correlation function of the fluctuating flux time series for each 10 min period. Then the time lag was determined where the auto-correlation has decreased by a factor of e−1 . Fig. 3 shows such an auto-correlation function for a period with wind speed of 5 m s−1 and a sensible heat flux of 280 W m−2 . Also shown are the auto-correlation function for the vertical wind speed (ρww ) which has a somewhat larger integral time scale, and an exponential function with decay time of 1.5 s. It is observed that the shape of ρ ff deviates from the assumed exponential form. Table 1 shows typical values for the integral time scales of w, the vertical velocity and for f, the flux fluctuation for both night Table 1 Integral scales for the vertical wind speed w and the flux fluctuation signal f for daytime and nighttimea
Fig. 2. Time series of 10 min averaged sensible heat fluxes based on sampling intervals 0.1 (full line), 1.0 (filled circles) and 4.0 s (plus) for 1 day.
Daytime Nighttime
Jw (s)
Jf (s)
Jw U/ze
Jf U/ze
4.0 1.8
2.5 1.3
1.0 0.3
0.6 0.2
a First two columns, time scales; last two columns, dimensionless scales.
F.C. Bosveld, A.C.M. Beljaars / Agricultural and Forest Meteorology 109 (2001) 39–45
43
Table 2 Nine days daytime (9:00–15:00 h) averaged sensible heat flux for the original time series (∆ = 0.1 s) and the stripped time seriesa Heat flux (W m−2 )
Average Bias
∆ (s) 0.1
0.2
0.4
1.0
202.64 –
202.62 −0.02 ± 0.05
202.96 0.33 ± 0.12
202.71 0.07 ± 0.30
2.0 200.80 −1.84 ± 0.59
4.0 202.12 −0.52 ± 1.15
a Bias is relative to the ∆ = 0.1 s time series. The standard error in the bias is derived from the modelled sample rate induced variance in the individual fluxes.
and daytime. It is observed that the integral time scale for the flux signal is considerably smaller than for the vertical wind speed. At nighttime integral time scales are smaller because the stable stratification limits the vertical movement of the air. Table 1 also shows the more fundamental dimensionless integral time scales. These values are obtained by scaling J with the wind speed U and the effective measuring height ze . In the absence of a continuous time series we cannot determine the natural variance, to assess the sample induced variances of the six time series. We can however calculate the sample induced variance of the five constructed fluctuating flux time series relative to the original 10 Hz time. Here we will show that, under the assumption T Jf , this variance is equal to the difference of the variances of the two individual time series, i.e. {(fjT (p∆) − fjT (∆))2 } = VfT (p∆) − VfT (∆)
(13)
where fjT (p∆) represents the time series with each pth sample retained. The proof follows from the observation that this is equivalent to {fjT (p∆)fjT (∆)} = {fjT (∆)fjT (∆)}
(14)
Thus, the ensemble covariance between the stripped time series and the original sampled time series has to be equal to the variance of the original sampled time series. For an individual sample of the original time series its covariance with this time series is equal to the variance of the original time series if the sample is sufficiently remote from the beginning and the end of the time series. Now, the stripped time series fjT (p∆) is nothing but a linear combination of such individual. This then proves Eq. (14) and hence Eq. (13). Edge effects become negligible when T Jf .
3.1. Finite sampling rate induced flux bias From the theoretical analysis, it follows that no bias is induced by a finite sampling interval. To test this hypothesis the averaged sensible heat flux was calculated over the 9 mid-day periods (9:00–15:00 UTC) for the different sampling intervals ∆. Table 2 gives the averaged values for each sampling interval. Also shown are the biases relative to the ∆ = 0.1 s case. For all cases biases are less than 1% of the averaged value. The hypothesis is that these biases are not significant. To test this hypothesis the standard error (68% confidence interval) in the biases were estimated from the modelled variances in the individual fluxes according to the Eqs. (12) and (13). It is observed that the biases are of the same order as the standard error, although for ∆ = 0.4 and 2.0 s the biases are quite large. These results give a reasonable confidence in the theoretical treatment of the bias. 3.2. Sampling rate induced variance in the flux From the theoretical analysis it follows that the variance in individual 10 min flux estimates increases with increasing sample interval. For each day and each of the five stripped time series, experimental values of the sampling rate induced flux variances were derived from the 144 values of FT (∆) − FT (0.1). Sampling rate induced model variances are derived by assuming an exponential auto-correlation function (Eqs. (11) and (12)). Fig. 4 compares the modelled values (open symbols) against observations. Standard deviations are shown calculated as the square root of the variances. For each sampling interval there are nine points in the figure corresponding to the nine different days. Values are displayed along logarithmic axis. At the top right side of the graph the total standard deviation
44
F.C. Bosveld, A.C.M. Beljaars / Agricultural and Forest Meteorology 109 (2001) 39–45
Fig. 4. Daily averaged sampling rate induced flux standard deviation for sampling intervals 0.2, 0.4, 1.0, 2.0 and 4.0 s minus the sampling rate induced flux standard deviation for sampling interval 0.1 s. Shown are model against observation values. Also shown is the natural standard deviation derived from the 0.1 s time series (observation) and from Eq. (12) (modelled). Open symbols represent model values based on an exponential auto-correlation function, filled symbols are based on the k = 0 term in Eq. (9).
for each day is shown for the ∆ = 0.1 s. case. These are derived from the sum of the natural variance and the sampling rate induced variance. By using the time series with the highest time resolution the latter variance is minimised. The observed total variance is derived by calculating the differences between subsequent 10 min fluxes. By taking 0.5 times the variance of the difference series, the variance in individual values is obtained, provided that the change in expected values for subsequent time intervals is small and provided that the statistical deviations are independent. The model values for the total variance are obtained from Eq. (12) by substituting the observed values for µf and Jf . The modelled total variances are somewhat smaller than the observed values. This might by caused by contributions to the observed natural variance on long time scales not captured by the assumed form of the auto-correlation function in the model calculations. The figure shows a reasonable agreement between modelled and observed sampling rate induced variances. In all cases modelled values underestimate the observed variances, although the relative differences seems to decrease for larger sampling intervals. As shown in Fig. 3, the model assumption of an exponen-
tial auto-correlation function is not entirely correct. The actual auto-correlation of the flux fluctuation signal is steeper at small delay times and less steep for delay time of the order Jf and larger. This will lead to an underestimation of the modelled variances for the shorter time intervals. A new calculation was performed by using the actual auto-correlation function for each 10 min period and calculating the sampling rate induced variance from the k = 0 term in Eq. (9) only (filled symbols in Fig. 4). In this case a better agreement between observed and modelled sampling rate induced variance was obtained. Fig. 4 also shows that the natural variance is larger than the sampling rate induced variance for all the sampling intervals used in this study. For a sampling rate of 1 s the sampling rate induced variance is on the average a factor 50 smaller than the natural variance.
4. Conclusions In this study, the impact of sampling on eddycovariance fluxes is investigated. It is shown both theoretically and experimentally that sampling does not introduce a bias in the fluxes. The analysis shows further that the current theory, with an assumed exponential auto-correlation function, can explain the sampling rate induced variance in the fluxes reasonably well. The sampling rate induced variance depends on the integral time scale of the flux fluctuation signal. This integral time scale depends on stability, the wind speed and the effective measuring height. For a sample interval of 1 s, an effective measuring height of 18 m and a typical daytime wind speed of 3.5 m s−1 , as described in Bosveld and Bouten (2001), it follows that the sample rate induced variance in the flux is approximately 50 times smaller than the natural variance in the flux. The current theoretical treatment shows that it is not the integral time scale of the auto-correlation function that determines the sampling rate induced variance in the flux, but the steepness of the auto-correlation functions between t = 0 and t = ∆/2. It is shown that the assumption of an exponential decaying auto-correlation function for the flux fluctuation signal is not correct. This leads to an underestimation of the sampling rate induced variance for sampling intervals smaller than the integral time scale.
F.C. Bosveld, A.C.M. Beljaars / Agricultural and Forest Meteorology 109 (2001) 39–45
References Bosveld, F.C., Bouten, W., 2001. Comparing transpiration models with eddy-correlation observations for a Douglas-fir forest. Agric. For. Meteorol. 108, 247–264. Haugen, D.A., 1978.I Effects of sampling rates and averaging periods on meteorological measurements. In: Proceedings of the Fourth Symposium on Meteorological Observations and Instrumentation, April 10–14, 1978, Denver, CO., American Meteorological Society, pp. 15–18. Lenschow, D.H., Mann, J., Kristensen, L., 1994. How long is long enough when measuring fluxes and other turbulence statistics? J. Atmos. Oceanic Technol. 11, 661–673.
45
Massman, W.J., 2000. A simple method for estimating frequency response corrections for eddy covariance systems. Agric. For. Meteorol. 104, 185–198. Moncrieff, J.B., Massheder, J.M., de Bruin, H., Elbers, J., Friborg, T., Heusinkveld, B., Kabat, P., Scott, S., Soegaard, H., Verhoef, A., 1997. A system to measure surface fluxes of momentum, sensible heat. J. Hydrol. 188/189, 589–611. Moore, C.J., 1986. Frequency response corrections for eddy correlation systems. Boundary Layer Meteorol. 37, 17–35. Wyngaard, J.C., 1973. On surface layer turbulence. In: Haugen, D.A. (Ed.), Proceeding of the Workshop on Micrometeorology, American Meteorological, Society, Science Press, Ephrata, PA, USA.
The impact of sampling rate on eddy-covariance flux estimates F.C. Bosveld a,∗ , A.C.M. Beljaars b a
Royal Netherlands Meteorological Institute, P.O. Box 201 3730 De Bilt, AE, The Netherlands b European Centre for Medium Range Weather Forecast, Reading, UK
Received 18 September 2000; received in revised form 9 May 2001; accepted 15 May 2001
Abstract It is shown that sampling rate has no influence on the expected value of an eddy-covariance flux estimate. The uncertainty in the flux estimate, however, increases with decreasing sampling rate. This sampling induced uncertainty is compared with the natural uncertainty in the flux estimate due to the stochastic nature of turbulence. It is shown that the sampling induced uncertainty becomes important when the sampling interval is larger than the integral time scale of the flux time series. Theoretical expressions are evaluated by using 9 days of turbulence data taken over a coniferous forest. For sampling intervals smaller than the integral time scale the sampling induced uncertainties are underestimated by the model because the actual auto-correlation function for the flux time series deviates from the exponential form assumed in the model. © 2001 Elsevier Science B.V. All rights reserved. Keywords: Scalar flux; Eddy-covariance; Integral time scale; Turbulence statistics
1. Introduction Vertical fluxes between the earth surface and the atmosphere can be measured directly in the atmospheric surface layer by means of the eddy-covariance method, provided that fast response sensors are available for the determination of the vertical wind speed and the transported quantity under consideration. Over the years significant improvements have been made in correction procedures because of limitations of the instruments and on the minimal requirements of instruments to give reliable flux estimates (Moore, 1986; ∗ Corresponding author. Tel.: +31-30-2206-911; fax: +31-30-2210-407. E-mail address: [email protected] (F.C. Bosveld).
Moncrieff et al., 1997; Massman, 2000). In these studies no attention is given to the subject of sampling rate requirements. The sampling rate required in an experiment depends on the purpose of the study. If one is interested in spectral analysis of the high frequency turbulent regime one is likely to choose the sampling rate related to the time constant of the measuring system. For example, when using a sonic anemometer with a path length of d at a wind speed U one is likely to use a sample rate faster than U/d. Another example is the shift of samples in time to account for delays in measuring tubes (i.e. trace gas analysers). One then would like to sample at a rate with which an accurate shift can be obtained. In such cases required sample rates are typically between 10 and 50 Hz. Note that for spectral purposes low-pass filtering at the Nyquist
0168-1923/01/$ – see front matter © 2001 Elsevier Science B.V. All rights reserved. PII: S 0 1 6 8 - 1 9 2 3 ( 0 1 ) 0 0 2 5 7 - X
40
F.C. Bosveld, A.C.M. Beljaars / Agricultural and Forest Meteorology 109 (2001) 39–45
frequency is desirable to avoid aliasing, whereas for flux estimates this filtering is undesirable due to loss of high frequency flux contribution. For monitoring of vertical fluxes one is often only interested in time averaged fluxes and much lower sampling rates may be chosen. Bosveld and Bouten (2001), for example, use a sampling rate of 1 Hz for eddy-covariance measurements taken over an 18 m high coniferous forest at a height of 30 m. Here we investigate the impact of sampling rate on time averaged flux estimates. Both the impact on the average flux and on the uncertainty of the flux estimate will be quantified. It seems that Haugen (1978) was the first to treat the subject of sampling requirements. Here we follow the treatment of Lenschow et al. (1994) who give a fairly complete analysis. We modify their analysis to avoid interference with the impact of the finite length of the time series, which they treat simultaneously. The latter phenomenon is related to uncertainties in the low-frequency contribution to the flux and is of no concern here. By treating sampling rate separately the nature of its impact on eddy-covariance fluxes will become more clear. The theoretical results are evaluated with turbulent measurements taken at the same location and height as those of Bosveld and Bouten (2001), but with a sampling rate of 10 Hz. Time series with lower sampling rates are constructed from the 10 Hz time series by retaining samples at regular intervals. Fluxes from this constructed time series are then compared with fluxes derived from the original sampled time series. In the following we will assume that spectral losses due to finite response characteristics of the measuring system can be ignored.
First, we derive the expected flux and show that an unbiased estimate is obtained independent of the sampling rate. For the cross product we have {wj (t)sj (t)} = F(t) = F
where F is the flux. The constantness of F follows from the stationarity of the time series over the period T. Each simultaneous measurement of w and s results in an unbiased estimate of the flux and thus the flux calculated from a sampled time series will result also in an unbiased estimate independent of the sampling rate. Note that in practical applications the average of wj (t) and sj (t) over the period T are subtracted from the time series. This results in an extra term in the expression for the flux, representing the loss in the co-variance due to the omission of time scales T and longer. It can be shown that this low frequency term is almost independent of sampling rate, provided that the sampling interval is much smaller than T. In this study, we are only interested in a comparison between time series with different sampling rates and we can omit this term in the analysis. Now, we turn to the relation between the sampling rate and the variance of the flux estimate. For a realisation j an estimate of the flux over the period T is obtained by 1 T FTj = wj (t)sj (t) dt (3) T 0 The expected value of FTj is of course the flux F. Let us define the flux deviation time series as fj (t) = wj (t)sj (t) − F
2. Theory Let wj (t) and sj (t) be time series for t = 0 to T of the vertical wind speed and a scalar quantity, respectively. The time series are supposed to be stationary stochastic processes drawn from an ensemble. The index j indicates the jth realisation of the ensemble. Furthermore, we assume that the expected values (defined as the ensemble average and represented by {·}) of w and s are 0, thus w and s represent the fluctuations. In formula {wj (t)} = 0,
{sj (t)} = 0
(1)
(2)
(4)
This is a stochastic function with an expected value 0 at each instant of time. We now construct a sampled time series with a sampling interval ∆. The averaging period T is divided in N = T /∆ intervals. Samples are taken halfway the sample intervals. This is not a principle point, but this form gives a more elegant transition to the integral formulation. The deviation from F over a period T for the jth realisation of the ensemble when sampled with a time interval ∆ is fjT (∆)
N−1 1 1 = fj n∆ + ∆ N 2 n=0
(5)
F.C. Bosveld, A.C.M. Beljaars / Agricultural and Forest Meteorology 109 (2001) 39–45
41
In the limit for ∆ approaching 0, this transforms into the integral over the period T of the continuous function. Here we are interested in the ensemble variance VfT (∆) of fjT (∆). To this end we write VfT (∆) = {(fjT (∆))2 } N−1 1 1 1 = 2 fj n∆+ ∆ fj m∆+ ∆ 2 2 N n,m=0
(6) This can be written in terms of the auto-correlation function ρff and the variance µf of the stochastic function f N−1 µf ρff ([n − m]∆) VfT (∆) = 2 N
(7)
n,m=0
and with k = n − m and T = N ∆, we obtain VfT (∆) =
µf T
N
ρff (k∆) 1 −
k=−N
|k∆| ∆ T
(8)
and thus, if T is much larger than the time scale on which ρff approaches 0 we can write VfT (∆) =
∞ µf ρff (k∆)∆ T
(9)
k=−∞
By letting ∆ approaching 0 we obtain Eq. (46) of Lenschow et al. (1994). µf ∞ 2µf Jf T Vf (0) = ρff (t) dt = (10) T −∞ T where Jf is the integral time scale by definition. We call this the natural variance of the flux, i.e. the variance induced by the stochastic nature of the turbulence. The sampling induced variance WfT (∆) is defined as the difference between the total variance (Eq. (9)) and the natural variance (Eq. (10)): WfT (∆) = VfT (∆) − VfT (0)
(11)
The summation of Eq. (9) is represented in Fig. 1 for an exponential auto-correlation function with integral time scale Jf and for three different sampling intervals, i.e. ∆ = Jf , Jf /2 and 0. The figure shows that the increase in variance of the sampled case
Fig. 1. The continuous (∆ = 0) and finite sample interval (∆ = Jf and Jf /2) auto-correlation function, ρff , as function of delay time scaled on the integral time scale (t/Jf ), representing the summation of Eq. (9).
compared to the continuous case is mainly caused by the k = 0 term. The other terms in the sampled case are good approximations for the integral of the continuous case in the respective intervals. In fact, it can be shown that for an exponential auto-correlation function the k = 0 terms contribute negatively to WfT (∆) for at most 13%. This maximum contribution is reached when ∆ is of the order of Jf . The magnitude of the k = 0 term decreases with decreasing sampling interval. The figure also shows that, it is not the integral time scale that determines sampling induced variance, but the behaviour of the correlation function at delay times between 0 and ∆/2. For an exponential auto-correlation function with integral time scale Jf an exact form is retrieved by Lenschow et al. (1994), Eq. (58) 2µf Jf ∆ ∆ VfT (∆) = (12) coth T 2Jf 2Jf From this equation we derive that for an exponential auto-correlation function and ∆ equal to Jf the sampling induced variance WfT (∆) is 10% of the natural variance.
3. Experimental evaluation To test the theory on the influence of sampling rate on the mean flux and the variance in the flux,
42
F.C. Bosveld, A.C.M. Beljaars / Agricultural and Forest Meteorology 109 (2001) 39–45
we used a series of 9 days of 10 Hz turbulence data of the vertical wind speed and the sonic temperature, both measured with a Kaijo Denki DAT-300 sonic anemometer–thermometer system. The data were obtained during August 1996 at the Speulderbos research site, a coniferous forest site near the village of Garderen, The Netherlands. The geographical co-ordinates are 52◦ 15 N, 5◦ 41 E and 52 m altitude. Site and instrumentation are described in Bosveld and Bouten (2001). Eddy-covariance measurements were performed at 31 m above the forest floor. With an average tree height at that time of 22 m and a thinned stand the displacement height is estimated at 14 m, resulting in an effective measuring height of 17 m. This value is 2 m smaller than the effective measuring height during the 1989 experiment described in Bosveld and Bouten (2001). Six time series with different sample intervals were obtained by omitting observations. The sample intervals created were 0.1, 0.2, 0.4, 1.0, 2.0 and 4.0 s. The data were processed over 10 min. Sensible heat flux was derived by calculating the covariance’s between the vertical wind speed and the temperature using the six time series with different sample intervals. Fig. 2 shows a time series of heat fluxes for one day for sampling intervals 0.1, 1.0, and 4.0 s. A clear increase in scatter is observed for the larger sampling intervals. To arrive at modelled flux variances from Eq. (12) the variance µf of the fluctuating flux time series and the integral time scale Jf were estimated for each 10 min interval of the 9 daytime series.
Fig. 3. Auto-correlation function of the flux fluctuating time series ρff and of the vertical wind speed ρ ww based on a time period with wind speed of 5 m s−1 and a sensible heat flux of 280 W m−2 . Also shown is an exponential decay function with time constant 1.5 s.
Wyngaard (1973), based on the Kanses experiment, related µf to the flux F. He found µf ≈ 4 F2 for daytime cases and µf ≈ 10 F2 for nighttime cases. These relations were tested for the current data set and confirmed. Together with Eq. (10) we derive for a wind speed of 3.5 m s−1 during daytime a typical uncertainty of 20% in a 10 min averaged flux. Jf was estimated by calculating the auto-correlation function of the fluctuating flux time series for each 10 min period. Then the time lag was determined where the auto-correlation has decreased by a factor of e−1 . Fig. 3 shows such an auto-correlation function for a period with wind speed of 5 m s−1 and a sensible heat flux of 280 W m−2 . Also shown are the auto-correlation function for the vertical wind speed (ρww ) which has a somewhat larger integral time scale, and an exponential function with decay time of 1.5 s. It is observed that the shape of ρ ff deviates from the assumed exponential form. Table 1 shows typical values for the integral time scales of w, the vertical velocity and for f, the flux fluctuation for both night Table 1 Integral scales for the vertical wind speed w and the flux fluctuation signal f for daytime and nighttimea
Fig. 2. Time series of 10 min averaged sensible heat fluxes based on sampling intervals 0.1 (full line), 1.0 (filled circles) and 4.0 s (plus) for 1 day.
Daytime Nighttime
Jw (s)
Jf (s)
Jw U/ze
Jf U/ze
4.0 1.8
2.5 1.3
1.0 0.3
0.6 0.2
a First two columns, time scales; last two columns, dimensionless scales.
F.C. Bosveld, A.C.M. Beljaars / Agricultural and Forest Meteorology 109 (2001) 39–45
43
Table 2 Nine days daytime (9:00–15:00 h) averaged sensible heat flux for the original time series (∆ = 0.1 s) and the stripped time seriesa Heat flux (W m−2 )
Average Bias
∆ (s) 0.1
0.2
0.4
1.0
202.64 –
202.62 −0.02 ± 0.05
202.96 0.33 ± 0.12
202.71 0.07 ± 0.30
2.0 200.80 −1.84 ± 0.59
4.0 202.12 −0.52 ± 1.15
a Bias is relative to the ∆ = 0.1 s time series. The standard error in the bias is derived from the modelled sample rate induced variance in the individual fluxes.
and daytime. It is observed that the integral time scale for the flux signal is considerably smaller than for the vertical wind speed. At nighttime integral time scales are smaller because the stable stratification limits the vertical movement of the air. Table 1 also shows the more fundamental dimensionless integral time scales. These values are obtained by scaling J with the wind speed U and the effective measuring height ze . In the absence of a continuous time series we cannot determine the natural variance, to assess the sample induced variances of the six time series. We can however calculate the sample induced variance of the five constructed fluctuating flux time series relative to the original 10 Hz time. Here we will show that, under the assumption T Jf , this variance is equal to the difference of the variances of the two individual time series, i.e. {(fjT (p∆) − fjT (∆))2 } = VfT (p∆) − VfT (∆)
(13)
where fjT (p∆) represents the time series with each pth sample retained. The proof follows from the observation that this is equivalent to {fjT (p∆)fjT (∆)} = {fjT (∆)fjT (∆)}
(14)
Thus, the ensemble covariance between the stripped time series and the original sampled time series has to be equal to the variance of the original sampled time series. For an individual sample of the original time series its covariance with this time series is equal to the variance of the original time series if the sample is sufficiently remote from the beginning and the end of the time series. Now, the stripped time series fjT (p∆) is nothing but a linear combination of such individual. This then proves Eq. (14) and hence Eq. (13). Edge effects become negligible when T Jf .
3.1. Finite sampling rate induced flux bias From the theoretical analysis, it follows that no bias is induced by a finite sampling interval. To test this hypothesis the averaged sensible heat flux was calculated over the 9 mid-day periods (9:00–15:00 UTC) for the different sampling intervals ∆. Table 2 gives the averaged values for each sampling interval. Also shown are the biases relative to the ∆ = 0.1 s case. For all cases biases are less than 1% of the averaged value. The hypothesis is that these biases are not significant. To test this hypothesis the standard error (68% confidence interval) in the biases were estimated from the modelled variances in the individual fluxes according to the Eqs. (12) and (13). It is observed that the biases are of the same order as the standard error, although for ∆ = 0.4 and 2.0 s the biases are quite large. These results give a reasonable confidence in the theoretical treatment of the bias. 3.2. Sampling rate induced variance in the flux From the theoretical analysis it follows that the variance in individual 10 min flux estimates increases with increasing sample interval. For each day and each of the five stripped time series, experimental values of the sampling rate induced flux variances were derived from the 144 values of FT (∆) − FT (0.1). Sampling rate induced model variances are derived by assuming an exponential auto-correlation function (Eqs. (11) and (12)). Fig. 4 compares the modelled values (open symbols) against observations. Standard deviations are shown calculated as the square root of the variances. For each sampling interval there are nine points in the figure corresponding to the nine different days. Values are displayed along logarithmic axis. At the top right side of the graph the total standard deviation
44
F.C. Bosveld, A.C.M. Beljaars / Agricultural and Forest Meteorology 109 (2001) 39–45
Fig. 4. Daily averaged sampling rate induced flux standard deviation for sampling intervals 0.2, 0.4, 1.0, 2.0 and 4.0 s minus the sampling rate induced flux standard deviation for sampling interval 0.1 s. Shown are model against observation values. Also shown is the natural standard deviation derived from the 0.1 s time series (observation) and from Eq. (12) (modelled). Open symbols represent model values based on an exponential auto-correlation function, filled symbols are based on the k = 0 term in Eq. (9).
for each day is shown for the ∆ = 0.1 s. case. These are derived from the sum of the natural variance and the sampling rate induced variance. By using the time series with the highest time resolution the latter variance is minimised. The observed total variance is derived by calculating the differences between subsequent 10 min fluxes. By taking 0.5 times the variance of the difference series, the variance in individual values is obtained, provided that the change in expected values for subsequent time intervals is small and provided that the statistical deviations are independent. The model values for the total variance are obtained from Eq. (12) by substituting the observed values for µf and Jf . The modelled total variances are somewhat smaller than the observed values. This might by caused by contributions to the observed natural variance on long time scales not captured by the assumed form of the auto-correlation function in the model calculations. The figure shows a reasonable agreement between modelled and observed sampling rate induced variances. In all cases modelled values underestimate the observed variances, although the relative differences seems to decrease for larger sampling intervals. As shown in Fig. 3, the model assumption of an exponen-
tial auto-correlation function is not entirely correct. The actual auto-correlation of the flux fluctuation signal is steeper at small delay times and less steep for delay time of the order Jf and larger. This will lead to an underestimation of the modelled variances for the shorter time intervals. A new calculation was performed by using the actual auto-correlation function for each 10 min period and calculating the sampling rate induced variance from the k = 0 term in Eq. (9) only (filled symbols in Fig. 4). In this case a better agreement between observed and modelled sampling rate induced variance was obtained. Fig. 4 also shows that the natural variance is larger than the sampling rate induced variance for all the sampling intervals used in this study. For a sampling rate of 1 s the sampling rate induced variance is on the average a factor 50 smaller than the natural variance.
4. Conclusions In this study, the impact of sampling on eddycovariance fluxes is investigated. It is shown both theoretically and experimentally that sampling does not introduce a bias in the fluxes. The analysis shows further that the current theory, with an assumed exponential auto-correlation function, can explain the sampling rate induced variance in the fluxes reasonably well. The sampling rate induced variance depends on the integral time scale of the flux fluctuation signal. This integral time scale depends on stability, the wind speed and the effective measuring height. For a sample interval of 1 s, an effective measuring height of 18 m and a typical daytime wind speed of 3.5 m s−1 , as described in Bosveld and Bouten (2001), it follows that the sample rate induced variance in the flux is approximately 50 times smaller than the natural variance in the flux. The current theoretical treatment shows that it is not the integral time scale of the auto-correlation function that determines the sampling rate induced variance in the flux, but the steepness of the auto-correlation functions between t = 0 and t = ∆/2. It is shown that the assumption of an exponential decaying auto-correlation function for the flux fluctuation signal is not correct. This leads to an underestimation of the sampling rate induced variance for sampling intervals smaller than the integral time scale.
F.C. Bosveld, A.C.M. Beljaars / Agricultural and Forest Meteorology 109 (2001) 39–45
References Bosveld, F.C., Bouten, W., 2001. Comparing transpiration models with eddy-correlation observations for a Douglas-fir forest. Agric. For. Meteorol. 108, 247–264. Haugen, D.A., 1978.I Effects of sampling rates and averaging periods on meteorological measurements. In: Proceedings of the Fourth Symposium on Meteorological Observations and Instrumentation, April 10–14, 1978, Denver, CO., American Meteorological Society, pp. 15–18. Lenschow, D.H., Mann, J., Kristensen, L., 1994. How long is long enough when measuring fluxes and other turbulence statistics? J. Atmos. Oceanic Technol. 11, 661–673.
45
Massman, W.J., 2000. A simple method for estimating frequency response corrections for eddy covariance systems. Agric. For. Meteorol. 104, 185–198. Moncrieff, J.B., Massheder, J.M., de Bruin, H., Elbers, J., Friborg, T., Heusinkveld, B., Kabat, P., Scott, S., Soegaard, H., Verhoef, A., 1997. A system to measure surface fluxes of momentum, sensible heat. J. Hydrol. 188/189, 589–611. Moore, C.J., 1986. Frequency response corrections for eddy correlation systems. Boundary Layer Meteorol. 37, 17–35. Wyngaard, J.C., 1973. On surface layer turbulence. In: Haugen, D.A. (Ed.), Proceeding of the Workshop on Micrometeorology, American Meteorological, Society, Science Press, Ephrata, PA, USA.