The Effect of Coordinate Rotation on the Eddy Covariance Flux

한국농림기상학회지 제 권 제 호 ,

9

2

(2007)

Korean Journal of Agricultural and Forest Meteorology, Vol. 9, No. 2, (2007), pp. 100~108

The Effect of Coordinate Rotation on the Eddy Covariance Flux Estimation in a Hilly KoFlux Forest Catchment Renmin Yuan1, Minseok Kang1, Sungbin Park1, Jinkyu Hong2*, Dongho Lee1, and Joon Kim1

1 Dept. of Atmospheric Sciences, Yonsei Univ., Seoul, 120-749, Korea Lab. for Environmental Physics, Dept. of Crop & Soil Sciences, University of Georgia, 1101 Experiment Street, Griffin, GA 30228, USA

2

(Received March 29, 2007; Accepted June 11, 2007)

경사진

KoFlux

산림유역에서 에디공분산 플럭스 산출에 미치는 좌표회전의 효과

렌민 유안 ·강민석 ·박성빈 ·홍진규 *·이동호 ·김 준 1

1

1

2

1

연세대학교 대기과학과, 조지아 대학교 작물 토양 과학과 (2007년 3월 29일 접수; 2007년 6월 11일 수락)

1

1

2

ABSTRACT The Gwangneung KoFlux supersite, located in a rugged mountain region, is characterized by a low wind speed due to a mountain-valley circulation and rolling terrain. Therefore, it is essential to understand the effect of coordinate rotation on flux measurements by the eddy-covariance method. In this paper, we review the properties of three orthogonal coordinate frames (i.e., double, triple, and planar fit rotations) and apply to flux data observed at the Gwangneung supersite. The mean offset of vertical wind speed of sonic anemometer was inferred from the planar fit (PF) coordinate rotation, yielding the diurnal variation of about ± 0.05 ms−1. Double rotation ( = = 0 ) produced virtually the same turbulent fluxes of heat, water, and CO2 as those from the PF rotation under windy conditions. The former, however, resulted in large biases under calm conditions. The friction velocity, an important scaling parameter in the atmospheric surface layer, was more sensitive to the choice of coordinate rotation method. : KoFlux, Gwangneung supersite, Double rotation, Planar fit, Friction velocity v

w

Key words

INTRODUCTION Long term measurements of surface-air exchanges by the eddy covariance depend on the accurate application of this technique to diverse conditions of the study site. Over tall canopies and in complex terrain, for example, surface fluxes may contain significant contributions at much lower frequencies (i.e., at periods much longer than the typical averaging period of 15 - 30 min) than expected from classical studies (Finnigan ., 2003). Therefore, rotating coordinates every period (e.g.,< 30 minutes) may lead to a systematic underestimation of the et al

Corresponding Author: Jinkyu Hong([email protected])

surface exchange. Typical coordinate rotation has been used for aligning the sensor perpendicular to the earth gravity over a flat terrain surface every averaging time (Tanner and Thurtell, 1969; Wesely, 1970). However, this natural wind system has theoretical limitations at non-ideal sites because it assumes a 1-dimensional approximation to the surface-layer mass/energy balance such that the sensor tilt is the only source of mean vertical motion, thereby forcing mean vertical flows to be zero (Lee ., 2004). Particularly, in a mountain region, mesoscale circulation (e.g., mountain-valley circulation) and et

al

101

Korean Journal of Agricultural and Forest Meteorology, Vol. 9, No. 2

synoptic air systems can generate the mean vertical motion and lateral Reynolds stress (Baldocchi et al., 2000). Furthermore, because of finite averaging period, in the natural wind system, the random sampling error of surface fluxes is relatively larger under light wind conditions (Wilczak et al., 2001). Over a complex terrain, the rotation to the natural wind system is contaminated not only by a sensor location bias, electronic offset, and flow distortion, but also by advective fluxes. Recently, several rotation methods were proposed to consider the terrain variation at non-ideal sites (e.g., Paw U et al., 2000; Wilczak et al., 2001). These new rotations calculate the mean streamline from an ensemble of observation data over weeks or longer. In particular, planar fit (PF) rotation by Wilczak et al. (2001) statistically provides the sonic anemometer offset, and therefore the instrumental offset can be eliminated in the flux calculation. Furthermore, the PF rotation can be used in assessing the 2- and 3-dimensional flow field like vertical advections. The KoFlux program (http://www.koflux.org) was launched in 2001 for understanding the surface-air exchanges of energy, water, and CO2 in key ecosystems of Monsoon Asia and a supersite was built in Gwangneung forest (Kim et al., 2006). Noticeably, most terrestrial ecosystems are located in heterogeneous and complex terrain in Korea. The Gwangneung supersite is also located in a forested mountainous landscape. The observation and numerical modeling at the Gwangneung supersite revealed that a mountain-valley circulation is dominant and that non-zero mean vertical motion is possible with low mean wind speed (≤ 2 m s−1) at 20 m above the canopy top (Hong and Kim, 2005). Thus, it is critical to investigate the effects of coordinate frames on the surface flux estimation in the KoFlux program. This paper is organized as follows. First, we briefly review the three coordinate rotation methods, followed by the site description. Then, we examine an appropriate averaging time period for planar fit rotation and compare the computed turbulent fluxes from different coordinate frames. Finally, summary and conclusions are provided.

II. MATERIALS AND METHODS 2.1. Coordinate rotations

The most commonly applied technique for determining the angles necessary to place the sonic anemometer into a streamwise coordinate system (i.e., natural wind

system) involves double rotations (DR) (Tanner and Thurtell, 1969; Wesely, 1970; Kaimal and Finnigan, 1994; Lee et al, 2004). By double rotations, the x-axis is aligned to the mean flow and the z-axis is perpendicular to the underlying surface as a right-hand system. The DR aligns the x-axis with the mean wind vector at the end of each averaging period. Therefore, the natural wind system allows us to calculate the online surface fluxes based on the assumption that all atmospheric quantities vary only in the vertical direction. The first rotation sets ν = 0 so that wind components after the first rotation are given by: u1 cos θ sin θ 0 um (1) = v1 – sin θ cos θ 1 vm 0 0 1 wm w1 where u, v, and w are longitudinal, lateral, and vertical wind components, respectively; θ = tan−1 ( vm ⁄ um ) ; and the subscripts m and 1 denote wind components before

and after the first rotation, respectively. The second rotation sets w = 0 so that the final velocities are then given by: u2 cos φ 0 sin φ u1 v2 = 0 1 0 v1 w2 – sin φ 0 cos φ w1

(2)

–1

where φ = tan ( w1 ⁄ u1 ) and subscript 2 denotes wind components after the second rotation. Implicitly, the DR assumes that the mean wind direction is the same as the direction of Reynolds stress. That is, vw = 0 after the alignment of x-axis into mean flow. The third rotation (TR) was proposed to align mean flow into mean Reynolds stress. In this third rotation, the new y- and z-axes are rotated around x-axis until the cross-stream stress becomes zero (i.e., vw = 0 ) and the matrix for the third set of rotation equations then becomes: u3 1 0 0 u2 v3 = 0 cos ψ sin ψ v2 (3) w3 w 0 – sin ψ cos ψ 2 –1 2 v 2 w 2 where ψ = tan --------------and subscript 3 denotes wind 2 2 v2 – w2

components after the triple rotation (Kaimal and Finnigan, 1994). Using the triple rotation, McMillan (1988) obtained an improvement of surface flux estimation at a sloped site. However, it is apparent that the TR is only

Yuan et al.: The Effect of Coordinate Rotation on the Eddy Covariance Flux...

valid at an ideal site with no advection and has several flaws at a non-ideal site. We cannot expect the same direction between mean flow and stress over a complex terrain. In particular, at non-ideal sites, additional rotation makes degradation of data quality and gives unrealistic estimations (Finnigan et al., 2003; Finnigan, 2004; Lee et al., 2004). Double and triple rotations are intuitive and relatively easy to apply. In homogeneous surface, the DR or TR allocates the mean vertical wind to the sensor tilt and serves in leveling the sonic anemometer. Thus, they permit only 1-dimensional conditions of vertical gradient. However, 2- and 3-dimensional flows such as thermal circulation can generate nonzero mean vertical wind speed (Finnigan, 1999). In this case, over-rotation results in the systematic bias in the surface flux computation. Furthermore, the random sampling error of surface fluxes is relatively larger in the DR and TR under light wind conditions because of finite averaging time (Wilczak et al., 2001). In general, scalar fluxes are not particularly sensitive to tilt bias but momentum fluxes significantly depend on the coordinate frame (Lee et al., 2004). This indicates that coordinate rotation can impact nighttime data filtering using friction velocity and, in turn, the estimation of annual summation of the net ecosystem production. The planar fit rotation for tilt correction was introduced to overcome the limitations of the DR and TR by Wilczak et al. (2001). The PF rotation is also a righthand orthogonal coordinate frame and uses an ensemble of individual 30-minute (or 15-minute) records. That is, based on the measured mean wind vector during the whole experimental period, a fitted plane is obtained using multiple-linear regression. Therefore, the PF rotation can make a stable coordinate frame and minimize the over-rotation due to the flow distortion by instruments and tower frame. The mean vertical velocity on this plane is zero, but individual data run has non-zero mean vertical wind. The PF rotation starts from multiple-linear regression: (4) wm = b0 + b1 um + b2 vm where subscript m denotes wind components measured by a sonic anemometer; and b0, b1, and b2 are regression coefficients estimated from linear regression. The coefficient b0 is the mean offset error because it is extremely difficult to ‘zero’ the transducers on a sonic

102

anemometer to eliminate mean wind speed biases. The coefficient b1 is the weighted contribution of u to the mean vertical velocity, and the coefficient b2 is the weighted contribution of v to the mean vertical velocity. Using these coefficients, we can construct the pitch (α) and roll (β) angles for the PF rotation: u1 um cos α sin α sin β – sin α cos β v1 = 0 vm ⋅ cos β sin β w1 sin α – cos α sin β cos α cos β wm – b0

(5)

–1

where α = sin ( –b1 ⁄ b + b + 1 ) and β = tan (b2). We should note that a sonic anemometer should not be moved or leveled when the PF rotation is applied because of multiple linear regression used in the PF rotation. Lastly, we rotate these intermediate winds and the stress tensor for each run so that x-axis is set along the mean wind (i.e., v = 0 ). 2 1

2 2

u2 cos γ sin γ 0 u1 v2 = – sin γ cos γ 0 ⋅ v1 w2 0 0 1 w1

–1

(6)

where γ = tan –1( v1 ⁄ u1 ) . More detailed information can be found in Wilczak et al. (2001) and Hong and Kim (2002).

2.2. Site and measurements

The measurements were made at two flux towers in the KoFlux Gwangneung experimental forest located northeast of Seoul, Korea (37°45'25.37''N, 127°9'11.62''E) (Fig. 1). The main tower is located near the headwater of the hilly deciduous forest catchment. The hillsides are dominated by slopes of 10 20o, with a maximum slope of 51o. The whole catchment area is based on weathered gneiss and schist with an eastern aspect. Typical soil depth is from 0.4 to 0.8 m and the soil texture is mainly sandy loam (Lim et al., 2003). The site is dominated by an old natural forest of Quercus sp. and Carpinus sp. (80-200 years old) and the average canopy height is about 18 m. The eddy-covariance system was installed at 20 m and 40 m on the Main Tower. The data have been recorded since late 2001 and were sampled at 10 Hz and averaged every 30 min. Fluxes of sensible heat, latent heat, and momentum were measured by a 3-dimensional sonic anemometer (CSAT3, Campbell Scientific, Inc.)

103


Fig. 1. Topography of the KoFlux Gwangneung forest.

and an open-path infrared gas analyzer (CS-7500, LiCor, Inc.). In this study, we used the data recorded by the EC system at the 40 m height from 1 Jan. to 31 Dec. 2006. For data processing, we used the data processing program developed in the KoFlux network (Hong and Kim, 2002). Here we do not apply air density correction by Webb . (1980) to focus only on the coordinate rotation effect. More details about the site can be found at the KoFlux website (http:// www.koflux.org). et al

Fig. 2. Varition of tilt angle with the wind direction.

III. RESULTS AND DISCUSSION 3.1. Wind distribution and tilt angles

Figs. 2 and 3 show the tilt angle variation with wind direction and its frequency distribution, respectively. Prevailing wind is southeasterly in daytime and northwesterly in nighttime, and vertical motion is related with this diurnal wind circulation. It indicates that a mountain-valley circulation mainly controls the wind circulation at this site rather than synoptic scale flows.


104

Fig. 3. The frequency distribution of the wind direction in 2006.

α)

Fig. 5. The variation of the pitch (

β)

and roll (

angles.

The subscripts ‘30’ and ‘all’ denote each 30° sector and all wind directions, respectively.

Fig. 4. Tilt planes for the 90 - 120° and 270 - 300° sectors.

digital elevation model (DEM) and satellite images. In the original paper by Wilczak . (2001), the PF rotation was applied regardless of wind direction. In this study, however, we apply the PF rotation not only with all wind directions but also with different wind directions (with a range of 30o each) to consider the effect of rolling topography around the flux tower (Fig. 4). Our analysis shows that the rotation angles of the coordinate frames obtained by using all wind directions depart considerably from those obtained by the PF rotation with every 30o. Such a significant difference is expected if we consider the complex topography around the tower, and again emphasizes that an appropriate PF rotation should be applied with different wind sectors. In case of southeasterly conditions (e.g., dayet al

We note that a simple sinusoidal function is appropriate to fit the apparent terrain slope against wind direction (Fig. 2). Such apparent terrain slope is expected if we consider that the slope and aspect following the axis of the mountain-valley circulation is nearly constant. Fig.4 shows the tilt planes obtained from the PF rotations for the two prevailing wind directions: (1) 90 120° (typically in daytime) and (b) 270 - 300° (generally in nighttime). These tilt planes are the mean stream planes, on which the mean vertical velocities are zero and the aspect slopes constructed by the PF rotation correctly represent the real topography analysis from

105


time), the difference was smaller. On the other hand, for the wind directions from SE to NW (e.g., nighttime), the tilt angle was relatively larger, and we speculate that the changes of slope and aspect across the forest catchment directly influence the tilt plane constructed by the PF rotation because typical flux footprints are > 1 km in nighttime. 3.2. Effect of the PF rotation on surface fluxes

We investigated how many half-hourly data are necessary for getting stable rotation angles from the PF rotation (Figs. 6 - 8). First, the instrument error of w (b0) from the PF rotation was about -0.05 m s−1 in daytime when southeasterly flows were dominant. However, this instrument error of w was around +0.1 m s−1 in nighttime. This indicates that the instrument error computed from the PF rotation is a function of the instrument electronic bias as well as the flow conditions. This has important physical implications in the 2and 3-dimensional flow characteristics, which impact the estimation of the net ecosystem production, for

The variations of the tilt angle and mean offset error for the 90-120° sector. Fig. 6.

The variations of the tilt angle and mean offset error for the 270 - 300° sector. Fig. 7.

The variations of the tilt angle and mean offset error for the 180 - 210° sector. Fig. 8.

example, in a non-ideal site such as the Gwangneung site (Lee and Hu, 2002). The positive (negative) nonzero offset observed in this study could results in overestimation (underestimation) of the net ecosystem production in nighttime (daytime). We note that the pitch and roll angles quickly converge to a constant value after 100 to 150 half hourly data. These correspond to about one to two weeks of data collected for each prevailing wind direction (i.e., 90 - 120o, 270 - 300o), likely reflecting the synoptic scale influences on the PF rotation (Figs. 6 and 7). We note that the roll angle for the 90 - 120o sector shows a gradual change up to 600 half-hourly data, indicating the possible influence of a seasonal change in phenology (e.g., leaf area index). For other sectors of less frequent wind directions, it would take much longer periods (i.e., weeks to months) (Fig. 8). Fig. 9 shows the comparison of the surface fluxes computed with the DR and PF rotation. A one-to-one relationship between the DR and PF rotation presented for scalar fluxes (i.e., heat, water vapor, and CO2) shows virtually no difference (with R2>0.93) between the two methods under windy condition ( U ≥ 1 m s–1 ; see the closed circles in Fig. 9). In case of the friction velocity (u*), two rotation methods also show a good agreement within 5 % but the coefficient of determination is relatively low (R2 = 0.88) due to many outliers. Under calm conditions ( U < 1 m s–1 ), the fluxes with the DR deviate from those with the PF rotation (see open circles in Fig. 9), indicating that the flux computation with the PF rotation is particularly important in nighttime. Despite the fact that DR and PF rotations did not show much difference in surface fluxes at the site, a caution must be exercised in interpreting individual


106

Comparison of surface fluxes computed from the DR and PF rotations under two different wind speeds: (1) U ≥ m s−1 (●) and (2) (○) U < 1 m s−1. H, LE, Fc and u* are the sensible heat, latent heat and CO2 fluxes, and friction velocity, respectively. Fig. 9.

Effects of the DR and PF rotations on the diurnal variations of H, LE, and u* on (1) 8 July (with variable winds) and (2) 5 August (in calm conditions) in 2006. Fig. 10.

107


surface flux measurement. For example, under windy conditions ( U ≥ 1 m s–1 in Fig. 10), the diurnal patterns and magnitudes of H, LE, and u* from the two rotations show no significant differences. Under calm conditions (August 5), however, both fluxes and u* from the PF rotation produced more plausible magnitudes and patterns while the results from the DR were unrealistic (e.g., LE = 790 W m and u* = 0.83 m s− between 1400 and 1430 local time when U = 0.1 m s− and Rn = 512 W m− ). -2

1

1

2

IV. SUMMARY AND CONCLUDING REMARKS In this paper, we briefly reviewed three coordinate rotation methods for the eddy-covariance technique: (1) double rotation ( v = w = 0 ), (2) triple rotation ( v = w = v′ w ′ = 0 ), and (3) planar fit rotation, and investigated the impact of these coordinate frames on the surface flux estimation at the Gwangneung KoFlux supersite located in a complex terrain. The double and triple rotations enable us to calculate the online surface fluxes. However, the DR and TR implicitly allow only the vertical gradient of wind and scalars and therefore have theoretical limits over a complex terrain where the cross wind shear, v′ w ′ , is typically non-zero. Furthermore, in a mountainous region, mesoscale flow like a mountain-valley circulation can generate non-zero mean vertical wind. The PF rotation can overcome several deficiencies in the DR and TR, and also has relatively small random sampling errors in light wind conditions. The PF rotation provides the offset of mean vertical wind speed so that we can attempt to compute the vertical advection using the residual vertical wind introduced by the PF rotation. Our analysis showed that more than two weeks of data is necessary to properly apply the PF rotation at the Gwangneung supersite but the length of data depended on wind directions. We found that the PF rotation should be applied with each 30° sector of wind direction to consider the effect of the rolling topography around the flux tower. The offset of the mean vertical wind speed was about | 0.05 | m s− with the opposite sign between daytime and nighttime. Turbulent fluxes can be biased due to the properly tile angle caused by such mean offset. Future work should focus on the dependency of the mean vertical wind bias from the PF rotation on the topography, canopy phenology, and atmospheric stability. 1

Under windy conditions, the DR and PF rotation produced virtually the same scalar fluxes, suggesting that using the online surface fluxes from the DR may be used to check the current status of the instrumentation in the fiel. However, we note that the DR produced a relatively poor correlation for friction velocity. Friction velocity is not only an important scaling variable in the atmospheric boundary layer but also a critical parameter for the nighttime data filtering (i.e., u* correction). In relation to such biases between the DR and PF rotation, further scrutiny on their relationship with vertical and horizontal advection is currently in progress.

ACKNOWLEDGMENTS This study is supported by a grant (Code: 1-8-2) from Sustainable Water Resources Research Center of 21 Century Frontier Research Program, the Eco-Technopia 21 Project from the Ministry of Environment, and the BK21 program from the Ministry of Education and Human Resource Management of Korea. We thank the two anonymous reviewers for their constructive comments and suggestions. st

적 요

험준한 산지에 위치한 광릉 KoFlux 슈퍼사이트는 기복이 진 지형과 산골 순환으로 인해 풍속이 약한 것이 특징이다. 따라서 에디 공분산 방법을 이용한 플 럭스 관측에 좌표회전이 미치는 효과를 이해하는 것이 필수적이다. 본 논문에서는 세 가지의 직교 좌표 체제 (이중, 삼중 및 planar fit 회전)의 특징을 살펴 보고, 이들을 광릉 슈퍼사이트에서 관측된 플럭스 자료에 적 용하였다. ‘Planar Fit(PF)’ 좌표회전으로부터 유추해 낸 초음파 풍속계의 연직 풍속의 평균 차감은 대략 ± 0.05 m s 이었다. 바람이 부는 경우, 이중 회전 ( v = w = 0 )에서 산출된 열, 수증기 및 CO 의 난류 플럭스는 PF 회전에서 나온 결과와 거의 같았다. 그러 나 바람이 잔잔한 경우에는 이중 회전이 큰 값의 편 의를 나타내었다. 대기 지표층의 중요한 스케일링 모수 인 마찰 속도는 좌표 회전 방법의 선택에 따라 더 민 감하게 반응하였다. −1

2

REFERENCE Baldocchi, D. D., J. J. Finnigan, K. Wilson, K. T. Paw U, E. Falge, 2000: On measuring net ecosystem exchange


over tall vegetation on complex terrain, Boundary Layer Meterorology , 257-292. Finnigan, J. J., 1999: A Comment on the Paper by Lee (1998): On Micrometeorological Observations of Surface-Air Exchange over Tall Vegetation, Agricultural and Forest Meteorology , 55-64. Finnigan, J. J., R. Clement, Y. Malhi, R. Leuning, and H. A. Cleugh, 2003: A re-evaluation of long-term flux measurement techniques-Part II: Averaging and coordinate rotation, Boundary-Layer Meteorology , 1-48. Finnigan, J. J., 2004: A re-evaluation of long-term flux measurement techniques-Part II: Coordinate systems, Boundary-Layer Meteorology , 1-41. Kaimal, J. J., and J. J. Finnigan, 1994: Atmospheric Boundary Layer Flows: Their Structure and Measurement, Oxford University Press, New York. Kim, J., D. H. Lee, J. K. Hong, S. K. Kang, S. J. Kim, S. K. Moon, J. H. Lim, Y. H. Son, J. S. Lee, S. H. Kim, N. C. Woo, K. H. Kim, B. Y. Lee, B. L. Lee, and S. Kim, 2007: HydroKorea and CarboKorea: Cross-scale studies of ecohydrology and biogeochemistry in a heterogeneous and complex forest catchment of Korea, Ecological Research , 881-889. Hong, J. and J. Kim, 2002: On processing raw data from micrometeorological field experiments, Korean Journal of Agricultural and Forest Meteorology., , 119-126. Hong, J., D. H. Lee, and J. Kim, 2005: Lessons from FIFE on scaling of surface fluxes at Gwangneung forest site, 96

97

107

113

21

4

Korean Journal of Agricultural and Forest Meteorology

, 4-14. (in Korean with English abstract)

7

108

Lee, X., and X. Hu, 2002: Forest-air fluxes of carbon, water and energy over non-flat terrain, Boundary-Layer Meteorology , 277-301. Lee, X., W. Massman, and B. Law, 2004: Handbook of Micrometeorology, Kluwer Academic Publishers, Dordrecht, The Netherlands. McMillen, R. T., 1988: An Eddy Correlation Technique with Extended Applicability to Non-Simple Terrain, Boundary-Layer Meteorology , 231-245. Paw U, K. T., D. D. Baldocchi, T. P. Meyers, and K. B. Wilson, 2000: Correction of eddy-covariance measurements incorporating both advective effects and density fluxes, Boundary-Layer Meteorology , pp. 487-511. Tanner, C. B., and G. W. Thurtell, 1969: Anemoclinometer measurements of Reynolds stress and heat transport in the atmospheric surface layer, Research and Development Tech. Report ECOM 66-G22-F to the US Army Electronics Command, Dept. Soil Science, Univ. of Wisconsin, Madison, WI. Webb, E. K., G. I. Pearman, and R. Leuning, 1980: Correction of flux measurements for density effects due to heat and water vapor transfer, Quarterly Journal of Royal Meteorology Society , 85-100. Wesely, M., 1970: Eddy correlation measurements in the atmospheric surface layer over agricultural crops. Ph. D. Dissertation. University of Wisconsin, Madison, Wisconsin, US. Wilczak, J. M., S. P. Oncley, and S. Stage, 2001: Sonic anemometer tilt correction algorithms, Boundary-Layer Meteorology , 127-150. 103

43

97

106

99