Tolerance of eddy covariance flux measurement

Eddy covariance (EC) flux measurement is the most-used technique for observation of the fluxes of sensible heat, latent heat, and carbon dioxide between the land surface and the atmospheric boundary layer. Despite the availability of plentiful EC data from numerous research projects, it is difficult to make meaningful comparisons of EC at different sites, to validate the models used, and to integrate observed data with models because the uncertainties of the method are inadequately defined. We developed a method to evaluate the uncertainties of the EC method without the need to consider individual site specifications and flux characteristics. We showed that the fractional error φ of EC (i.e., tolerance T) can be separated into random and illegitimate components. T can be used as a scale parameter for spatiotemporal stationarity, and can be defined as a rational function. We demonstrated a practical application of T analyses for two contrasting areas: a low-relief paddy field and an area of more complex land-forms where dramatic wind changes affect fluxes, and showed that T analysis provides an appropriate and effective method to determine the uncertainties in EC.


INTRODUCTION
The uncertainty generated by fluctuations in measurements of a quantity can be defined as the measurement error δ, which can be decomposed into random error δ r due to instrumental uncertainties or statistical fluctuations, systematic error δ s originating from defective equipment calibrations or discrepancies among different observers, and illegitimate error δ i arising from mistakes or measurement blunders (Bevington and Robinson, 2003). Commonly, δ r , δ s , and δ i , are called noise, bias, and outliers, respectively.
In tower-based eddy covariance (EC) flux measurement to monitor heat, H 2 O, and CO 2 exchanges above the land surface (Baldocchi, 2003), δ is apportioned into δ r and δ s under the assumption that all data contaminated by δ i are discarded during quality control and quality assurance (QC/ QA) procedures based on micrometeorological knowledge (Moncrieff et al., 1996;Loescher et al., 2006;Vickers et al., 2009). Whereas Vickers and Mahrt (1997) indicated that mesoscale variability and inhomogeneity (i.e., nonstationarity) are also constituents of δ, along with δ r and δ s , and Mahrt (1998) suggested that nonstationarity should be specified in EC.
We hypothesize that nonstationarity is a component of δ i and is inherent in EC because δ i arises from the mismatches between actual circumstances and potential preconditions based on micrometeorological measurement theories, and then it is not noise or bias but the outlier generated by applicational mistake or blunder behind the supportable micrometeorological background. In other words, it is difficult to presume the perfect observational theory satisfying every measurement condition without precondition over micrometeorological knowledge ever since. Therefore, the aims of this study were (1) to demonstrate the above hypothesis in terms of the fractional error φ of EC (tolerance T hereafter) based on the work of Finkelstein and Sims (2001) and Kim et al. (2008Kim et al. ( , 2009 (2) to separate the fractional random error φ r and the fractional illegitimate error φ i under the assumption that instrumental bias has been corrected; and (3) to test the validity of this approach as a QC/QA tool for field EC.

Field data
For T analysis of a homogeneous site, we used a week EC measurement from a paddy field in Sukhothai, Thailand (PST: 17°03'51''N, 99°42'17''E, 50 m asl). The data were recorded with a three-dimensional sonic anemometer (CSAT3; Campbell Scientific, Utah, USA) and an open-path CO 2 /H 2 O gas analyzer (LI7500; LI-COR, Nebraska, USA); both were deployed 7 m above the ground. The topography at PST is flat with sufficient fetches (> 700 m) in all wind direction for experimental period. Two different growing stages of the paddy were cultivated separately, and those leaf area index (LAI) is 0.5 at the side of clockwise from −45 to +135° from the direction of the sonic head and 2.0 at the other side, respectively.
For T analysis of a heterogeneous site, we used a week EC measurement from an area of mixed land cover in Tak, Thailand (DTT: 16°56'24'7N, 99°25'48''E, 110 m asl), with the same instrumentation as that used at PST but operated at a height of 30 m. The fetch area at DTT is relatively flat Correspondence to: Wonsik Kim, National Institute for Agro-Environmental Sciences, 3-1-3 Kannondai, Tsukuba, Ibaraki 305-8604, Japan. E-mail: wonsik@affrc.go.jp ©2011, Japan Society of Hydrology and Water Resources.
with gently undulating hills. When the data were recorded, the land cover was a mosaic of 60% crop fields (mainly paddy, cassava, corn and non-arable fields) with LAI from 0 to about 3, and 40% monsoon forest dominated by around 10-m-tall Shorea with LAI around 2.

Tolerance
The T of EC measurement is defined as , ( where F is flux. If δ s is zero (i.e., no bias), and δ r has a linear relationship to |F| with slope φ r (i.e., δ r = φ r |F|), equation (1) can be modified to .
( 2 ) It is a rational function having |F| as an independent variable under assumptions that φ r is a constant and δ i is outlier in the function. Computational sequence to estimate T for a period by EC measurement is: First, hourly F and T is calculated by (Leuning, 2004, Equations 6.24, 6.26 and 6.27) and (Kim et al., 2009, Equation 1), respectively. In particular, Finkelstein and Sims (2001) is applicable to estimate T because of merely employing auto-and cross-covariance instead of any preconditional assumptions of time series, boundary layer states, and an autocorrelation function of integral time scale; second, expected T for a period according to the kind of F is estimated by mode which denotes the T point having most high frequency of hourly T among 0.01 T interval (frequency panels in Figure 1); third, using expected T as a coefficient and F as independent variables, a robust nonlinear fitting [R function nlrob (http://www.rproject.org/; http://svitsrv25.epfl.ch/R-doc/library/robustbase/ html/nlrob.html)] is assessed to determine δ i and scale parameter σ as one of deviation parameters in equation (2); finally, the minimum T value that produced the smallest σ is assigned to φ r in equation (2)

Tolerance
Estimated hourly T at PST and DTT fitted equation (2) T well as a function of F (Figure 1). The fit was better when weekly T was smaller and kurtosis was higher for hourly T (Figure 1). Therefore, T, as expressed in equation (2), is a rational function for the observational data we used. That is to say that consideration of only δ r in T analysis is inappropriate, because δ i is inherent in almost all EC measurement (e.g., Figure 3 in Finkelstein and Sims (2001), Figure 2 in Hollinger and Richardson (2005) and Figure 9 in Vickers et al. (2009)).
As outliers, the hourly T including δ i originated from temporal nonstationarity, principally appears where the values of F are near zero on almost all F types and at both sites, and the T having δ i due to spatial nonstationarity, predominantly reveals where the values of F are around lE > 50 W m −2 and F c < −0.1 mg m −2 s −1 at DTT (Figure 1). These locational differences are helpful for understanding whether variation of atmospheric background or heterogeneity of surface conditions cause the distortion of stationarity in those larger T.
Furthermore, considering that T was affected by heterogeneity of not only micrometeorological events (i.e., convection and horizontal divergence; Loescher et al., 2006) but also biological responses (i.e., photosynthesis and transpiration; Oren et al., 2006), and minimum T (i.e., φ r ≈ 0.07) could be a constant (Figure 1), therefore, we suggest that T is a general scale parameter of nonstationarity for EC measurement.

Fractional random error
If δ i is eliminated, theoretically from equation (2), T equals φ r and is constant. Alternatively, if EC measurement is observed under ideal experimental conditions, T represents φ r . Hence, the minimum value of T with the smallest value of σ in a periodic analysis changing the time span for the calculation of H, lE, and F c (Table I) can appropriately be defined as φ r for observational EC data. Our analysis showed that φ r was 0.07 ± 0.01. It could be the measurement value of φ r under optimized condition of EC because it is the smallest value compared with those of the other periods and studies.
Approximations such as that above have been reported in many previous studies (Lenschow et al., 1994;Finkelstein and Sims, 2001;Hollinger and Richardson, 2005;Richardson et al., 2008;Vickers et al., 2009); however, those studies did not provide enough persuasive data to reach the conclusion we have presented here because of a lack of clarity about how to classify the inherent nonstationarity in uncertainty analysis, and the inability to provide a general method to quantify it.
However, our successful determination of φ r (i.e., 0.07 ± 0.01, T lE of PST for 24-26 May 2007 in Table I) might account for the almost perfect spatial homogeneity by a specific cultivation subject to water body. Even at DTT of heterogeneous land cover, not only T Fc (0.09 ± 0.02 for 23-25 April 2009 in Table I) approached φ r value, but also some hourly T reached the determinate φ r value under the condition of friction velocity u * > 0.5 m s −1 (Figure 1). These results are of brace to understand our φ r determination in observational base.

Quality control and quality assurance
The T ± σ used as a tolerance criterion provides assured quality control of EC data. For example, for lE data acquired over one week at PST, 55% (circles in left middle flux panel of Figure 2) is acceptable with a tolerance level of 8% based on 67.8% probability; 34% (crosses in the same panel) are rejected; and 11% (shown in gray regions of the panel, not symbolized) are filtered before this analysis because two more consecutive spikes are detected over the turbulence trend. Even in the worst case, at DTT, 40% of F c data is acceptable at a tolerance level of 14%, 36% are rejected, and 24% are discarded (right lower flux panel in Figure 2). Conveniently, it is possible to estimate QC/QA criteria from T analysis for specified periods, sites and flux types as shown in Table I.
T analysis provides criteria not only for uncertainty, but also provides important QC/QA that is independent of site characteristics. That is, it overcomes the need for sitespecific investigations, such as footprint analysis and u * correction, because it includes uncertainty of spatiotemporal nonstationarity without precondition. For example, the estimated T included the fractional illegitimate error φ i attributable to distortion of stationarity at PST (1-3%) and DTT (2-9%) based on δ i ≈ (T − φ r )|F|, where φ r = 0.07 (see Fractional random error section).
Therefore, EC with uncertainty information based on T analysis can be used for comparison of sites, validation of models, and advanced forward investigations. To reserve higher quality than estimated uncertainty value of each flux, we recommend that the rejected measurements (crosses in Figure 2) are re-estimated by using a gap-filling method, or are screened by more strict critical value of T irrespective of data continuity. Furthermore, we found that the relationship between T and u * was not clear, thus it was difficult to say that every F value having u * > 0.2 m 2 s −1 always located in acceptable T estimated in this experimental period (see the gray regions in Figure 1).

CONCLUSIONS
We proposed a rational function T that can be used to determine fractional error φ of EC and demonstrated its validity for use with observational EC data. The function clearly shows that, T is composed of φ r , which is ever present in the randomness of white noise, and φ i , which varies according to the distortion of spatiotemporal stationarity of EC measurement. We demonstrated φ i to be constant with a value 0.07 even though additional studies are advisable, and successfully evaluated φ i from the discrepancy between T and φ r . Inferring from previous studies (Finkelstein and Sims, 2001;Hollinger and Richardson, 2005;Mano et al., 2007;Kim et al., 2011) it is therefore probable consequence that φ r by EC measurement under optimized micrometeorological condition is a constant because H, lE and F c are all controlled by eddy mixing process.
The T analysis is also useful for QC/QA of EC data. It provides information about uncertainty, a criterion for data filtering, and a scale parameter for nonstationarity, without the need for other analysis or correction. In addition to dealing with uncertainty and stationarity, T analysis is independent of site specific and flux kind.
For future investigations, we will first consider whether or not T values should be smaller than 7%. We will then attempt to develop a heterogeneity index using the considered T and φ i and, finally, develop a model that will help to reduce the uncertainty from future EC data analyses. Figure 2. Time series of estimated hourly fluxes of sensible heat H, latent heat lE and CO 2 F c showing tolerance T (whiskers), and associated micrometeorological parameters wind direction W d ; surface-layer scaling parameter zL −1 where z is measurement height and L denotes the scale of the Monin-Obukhov length in boundary layer; and friction velocity u * . Hourly F is shown by circles for weekly data within T ± σ (left-upper numeral) or by crosses for weekly data that are not within T ± σ. Data of acceptable F (circles) and the associated parameters are shaded according to W d . Vertical gray regions are the screen for F from consecutive spikes more than twice in its one-hour turbulence trend.