Local flux-profile relationships of wind speed and temperature in a canopy layer in atmospheric stable conditions

Introduction Conclusions References


Introduction
Fluxes of momentum and heat in the atmospheric surface layer are generally described with non-dimensional wind and temperature gradients using Monin-Obukhov similarity theory.These are usually expressed as functions of the atmospheric stability parameter z/L, where z is the measurement height and L is the Obukhov length (Webb, 1970;Businger et al., 1971;Dyer, 1974;Högström, 1996).Similarity functions agree with direct measurements for a wide range of unstable and neutral conditions (Högström, 1988;Sugita and Brutsaert, 1992).However, the agreement is not nearly as evident in moderately and strongly stable conditions (Lee, 1997;Sharan et al., 2003).
In stably stratified boundary layers, turbulence is generated by mechanical shear.It is destroyed by negative buoyancy forces and viscous dissipation, with the suppression of turbulence by negative buoyancy limiting the extent of vertical mixing (Stull, 1988).Nieuwstadt (1984a, b;1985) suggested that local Obukhov length ( ) based on local turbulent fluxes should be a more fundamental length scale in stable stratification.As a result, dimensionless gradients of wind speed and temperature could be expressed as functions of a single scaling parameter, or local stability parameter, ζ (= z/ ).Local scaling has been applied in studies of fluxprofile relationships in the stable atmospheric boundary layer (Howell and Sun, 1999;Grachev et al., 2005;Yagüe et al., 2006;Basu et al., 2006;van de Wiel et al., 2008).From direct numerical simulation of the Navier-Stokes equations, van de Wiel et al. (2008) showed that in true homogeneous and stationary conditions the similarity functions remain valid even in strong stability.On the other hand, it has been shown that, Published by Copernicus Publications on behalf of the European Geosciences Union.
G. Zhang et al.: Local flux-profile relationships of wind speed and temperature as the stability parameter ζ increases, the similarity functions tend to level off, i.e., in very stable conditions, stability tends not to control momentum and heat fluxes (Grachev et al., 2005;Yagüe et al., 2006).As indicated by Mahrt (2007), the "leveling off" is mainly due to non-stationarity and nonhomogeneity.
A forest canopy acts as a source and/or a sink of energy, carbon, and water vapor.An understanding of the momentum and heat transfer inside the canopy is essential to quantify and model the atmosphere-biosphere turbulent exchange.Simpson et al. (1998) showed that measured eddy-covariance fluxes agree well with those derived using similarity theory at about twice canopy height.They also noticed significant discrepancies at heights closer to the canopy (1.2 to 1.4 h, where h is the canopy height).These suggest that the Monin-Obukhov similarity theory work well above the roughness sublayer but not really as well within the roughness sublayer.
Other parameters, such as h/L (the ratio of the canopy height h to the Obukhov length L at the top of or above the canopy), have also been used to describe the influence of the atmospheric stability on turbulence properties inside the canopy (Shaw et al., 1988;Leclerc et al., 1990Leclerc et al., , 1991;;Brunet and Irvine, 2000).Using essentially the same approach as Shaw et al. (1988) and Leclerc et al. (1990Leclerc et al. ( , 1991)), Kruijt et al. (2000) analyzed within-canopy turbulence properties using five stability classes with 1/L.More recently, Harman and Finnigan (2007) defined a canopy penetration depth Lc and used Lc/L as the atmospheric stability parameter.However, thermal stratification within a canopy generally differs from that above the canopy.This is especially true in stable atmospheric conditions, as the flow is impacted by the canopy structure which modulates canopy turbulence and thermal structures.
Different stability parameters have been used to represent the stability within the canopy.Lee and Mahrt (2005) formulated a revised local stability parameter, γ z/ , to study the effect of stability on mixing in forest canopies, where is the local Obukhov length and the coefficient γ = (z r − d)/z r , with z r for the top of the roughness sublayer and d for the displacement height.Jacobs et al. (1992) showed that the bulk Richardson number within a maize crop canopy provides more useful information in calm nocturnal conditions than the Obukhov length either within or above the canopy.Leuning (2000) used h/L as stability parameter, but the Obukhov length L of the roughness sublayer at night was estimated from the above-canopy gradient Richardson number, not from the sensible heat flux and friction velocity measured with the eddy-covariance method.In this way, the fluxes estimated from the inverse Lagrangian analysis agreed well with eddy-covariance measured fluxes.This study suggests that, in stable conditions, the gradient Richardson number should be a better parameter than z/L when L is obtained directly from eddy-covariance measurements.Baas et al. (2006) also indicated that in the stable boundary layer, fluxes are generally small and gradients of mean variables large.The authors thus suggest that a gradient-based scaling might be more suitable than a flux-based scaling.The latter is usually used for the unstable boundary layers with large fluxes and small gradients.
In this paper, we study flux-profile relationships at different depths within a tall forest canopy, using eddy-covariance measurements and temperature profile data within and above the canopy.We first examine the profiles of wind speed and temperature, vertical gradients, and momentum and heat fluxes, all quantities involved in stability parameters and flux-profile relationships of wind and temperature.Then we compare the relationships of local non-dimensional wind and temperature gradients (φ m and φ h ) with different atmospheric stability parameters, i.e., h/L(h) at tree top, mean local z/L(z), and local bulk Richardson number Ri, at different layers within the forest canopy in nighttime stable atmospheric conditions indicated with h/L(h) at tree top.Both the non-dimensional gradients and stability parameter share one or more variables and self-correlation is not avoidable.Therefore, we use a self-correlation analysis (Mahrt and Vickers, 2003;Klipp and Mahrt, 2004;Ha et al., 2007;Mahrt, 2007) to check the degree of physical meaning involved in the above relations.

Experimental site
The experiment was conducted in a 15-year-old managed slash pine forest canopy (Pinus elliottii L.) at the Florida AmeriFlux site (29 • 45 N, 82 • 10 W) near Gainesville, FL.The site is flat and has an elevation of about 50 m.The reader is referred to Gholz and Clark (2002) for more climatic information about the site.
There were two peaks in the canopy leaf area index (LAI) density profile (Fig. 1, measured in June 2003, with canopy height of 13.1 m and total LAI of 2.85 ± 0.99), one at 0.75-0.8h in the tree crown and the second at 0.1 h corresponding to the understory.The understory was uneven and composed mainly of saw palmetto, wax myrtle, Carolina jasmine, gall berry, and eastern bracken fern.

Instrumentation and data collection
During the experimental campaign from 15 April to 11 July 2004, seven three-dimensional sonic anemometers (CSAT3, Campbell Scientific Inc., Logan, UT) were placed at 0.08, 0.21, 0.48, 0.59, 0.89, 1.0, 1.4 h, where the canopy height h was 13.5 m.The data were logged into a CR9000 data logger (Campbell Scientific Inc., Logan, UT) at 10 Hz.Fourteen aspirated thermocouple sensors (Campbell Scientific Inc., Logan, UT) were mounted on a nearby tower at levels 0.07, 0.15, 0.23, 0.3, 0.39, 0.45, 0.53, 0.6, 0.68, 0.76, 0.84, 0.91, 0.99, and 1.06 h to get a temperature profile.Halfhourly mean temperature was recorded with a CR-23X data logger (Campbell Scientific Inc., Logan, UT).In addition, an automatic weather station was installed above the forest to record half-hourly mean air temperature, wind speed and direction, relative humidity, solar and net radiance, atmospheric pressure and rainfall.

Data processing
Data were processed in 30-minute intervals.The spikes in the data were detected and removed with a method proposed by Højstrup (1993) and Vickers and Mahrt (1997).The planar fit (PF) method (Wilczak et al., 2001) was applied to correct sonic anemometer tilt errors.In addition, the temporal trend in each individual run was removed with linear regression before various statistical properties were calculated.
The local bulk Richardson number, Ri, was calculated using the potential temperature (estimated from the measured air temperature profile) and wind speed at two heights (z 1 and z 2 ) in a layer within the canopy: where g is the acceleration due to gravity, θ is 30-min mean potential temperature with θ representing the spatial mean potential temperature in the layer, and U is the 30-min mean cup wind speed (U = (u 2 + v 2 ) 1/2 ) calculated from the sonic anemometer wind components u and v.
In addition to the stability parameter h/L at tree top, the local stability parameter (z/L) was also calculated from the local Obukhov length L.
where k is the von Karman constant (taken as 0.4), θ v is the virtual potential temperature from sonic anemometer, w θ v is the local covariance between the vertical wind velocity and the virtual potential temperature, and u * is the local friction velocity.The local friction velocity u * (z) was calculated as Local non-dimensional wind and temperature gradients (φ m and φ h ) were calculated using gradients between the sonic anemometer levels within the canopy layer: and where θ * is the local temperature scale estimated using the local friction velocity and local covariance of the vertical velocity and temperature as The degree of self-correlation between the nondimensional gradients and the stability parameter due to sharing one or more variables was also estimated.We determined the self-correlation between φ m and Ri and between φ h and Ri using the method proposed by Mahrt and Vickers (2003) and by Klipp and Mahrt (2004).The procedures are: 1. each variable involved in the calculation of φ m , φ h and Ri is picked at random from the original data set; 2. φ m , φ h , and Ri are calculated with the randomly picked variables; 3. repeat the above-stated two steps as many times as the number of the original data; 4. in contrast with Mahrt and Vickers (2003) and Klipp and Mahrt (2004) where they used only a linear correlation, a linear correlation coefficient or the coefficient of determination (i.e., ratio of the variance of simulated data to the variance of data to be simulated) is calculated between φ m (or φ h ) and Ri, depending on either linear or nonlinear relationship between them in different ranges of Ri; and 5. the above steps are then repeated 1000 times and the coefficients averaged.
In this method, the measured values of variables contained in φ m , φ h and Ri, i.e., friction velocity, temperature scale, potential temperature, wind shear, and potential temperature gradient, are redistributed randomly.The φ m (or φ h ) and Ri created from the randomly redistributed variables will not retain any real connection with physical meaning.Thus, the correlation between them is a measure of self-correlation due to sharing one or several variables.The difference between the coefficient of determination of the original data set and that of the randomly created data set is an estimate of the physical variance explained by the examined relationship (Mahrt and Vickers, 2003;Klipp and Mahrt, 2004;Baas et al., 2006;Vickers et al., 2009).The present analysis is focused exclusively to the nighttime dataset (19:30-05:30 h EST) with stable conditions indicated by h/L(h) > 0.

Vertical profiles within the canopy
Vertical profiles of atmospheric variables used in calculating local stability parameters and non-dimensional wind and temperature gradients within the forest are plotted in Fig. 2. The data are grouped and averaged according to the atmospheric stability at tree top, h/L(h).Three classes are shown to represent a range of weakly, moderately and strongly stable conditions at nighttime.Table 1 shows the ranges of h/L(h) values of these classes along with their mean, standard deviation, and number of periods in each class.It should be pointed out that, for those variables at the level 0.21 h or in the layers 0.48-0.21h and 0.21-0.08h, the period numbers are fewer than those at upper levels due to the removal of occasional erroneous sonic anemometer data points at 0.21 h.

Wind speed gradient and momentum flux
The profiles of horizontal wind speed and its vertical gradient in different stability conditions are shown in Fig. 2a and b.As expected, the horizontal wind speed decreases with depth in the canopy and with stability.It is noticed that, in the moderately and strongly stable conditions, a local maximum of wind speed exists at about 0.5 h (Fig. 2a).This coincides with the region of the sub-canopy characterized by a sparse foliage just beneath the dense crown foliage.Wind gradients are the largest in the crown as 0.23±0.057,0.15±0.037,and 0.12 ± 0.059 s −1 for weakly, moderately, and strongly stable classes, respectively, and then rapidly decrease with depth until the mid-canopy layer (Fig. 2b).This is due to the absorption of momentum by leaves and branches, leading to a fast decrease of wind speed in the upper canopy.The gradients are small in lower layers of the canopy, with mean values as 0.022, 0.016, and 0.014 s −1 for weakly, moderately, and strongly stable classes, respectively.The wind speed gradient inside the canopy decreases with increasing stability, same as above the forest (data not shown).Negative velocity gradients are usually observed in moderately and strongly stable atmospheric conditions at about mid-canopy and corresponds to a local maximum of wind speed in the sparse layer of the trunk space.
Figure 2c presents the profiles of the non-dimensional covariance u w normalized with the friction velocity at tree top in different stability conditions.In the crown layer, the mean normalized covariance values change from about −1 at tree top to about zero at 0.6 h (about the bottom of the crown), indicating that the downward momentum flux decreases quickly with depth, due to the absorption by the crown.The mean normalized covariance is small and mostly positive (around 0.2 in mean) in the lower canopy in strongly stable conditions.This reflects the occurrence of upward momentum transport in strongly stable atmospheric conditions.The profiles of the normalized u w are similar to those of Kruijt et al. (2000).

Potential temperature gradient and heat flux
The profiles of potential temperature (Fig. 2d) and its vertical gradients (Fig. 2e) exhibit different patterns in different stable conditions.In weakly stable conditions, the potential temperature slowly decreases with depth in the upper canopy and decreases faster in the lower canopy (Fig. 2d).The gradients are small (around 0.02 Km −1 ) and decreases slowly with depth in the upper canopy.However, they rapidly increase with depth in the lower canopy to reach close to 0.4 Km −1 in mean value (Fig. 2e).The results in weakly stable conditions suggest an intrusion and mixing of air from above the canopy into the forest to about half of the canopy, with a temperature inversion deep inside the canopy due to the ground energy loss through long wave radiation.
However, in strongly stable conditions, the potential temperature decreases faster in the upper canopy region and slower in the forest floor (Fig. 2d).The temperature gradients are the largest in the top layer with a mean of 0.32 Km −1 and decrease with depth.They are close to zero at about the middle layer and then increase with depth until the lowest layer where they are 0.13 Km −1 on average (Fig. 2e).This suggests that the crown functions as a lid and limits the exchange of above-crown warmer air and below-crown colder air.There are thus several simultaneous distinct thermal stratification regimes in the crown, trunk, and understory layers.
The profiles of covariance of vertical velocity w and potential temperature θ, normalized with wind and temperature scales u * and |θ * | at tree top, are plotted in Fig. 2f, where the absolute value of θ * is used to maintain the sign of the local covariance w θ .In stable conditions, the normalized covariance is negative in the upper layer of the crown, which indicates downward heat transfer, the same direction as the temperature gradients in Fig. 2e.On average, it changes rapidly from −1 at the tree top to −0.12, 0.043, and 0.51 at the bottom of the crown for the weakly, moderately, and strongly stable atmospheric conditions, reflecting the rapid decrease in downward heat transfer in the crown.At mid-canopy and deeper inside the canopy, the mean normalized covariance changes little with height for weakly and moderately stable cases, with mean values of −0.12 to −0.08 and around zero, respectively.In strongly stable conditions, however, it decreases from 0.7 at 0.48 h to 0.3 at 0.08 h.By comparing Fig. 2e and f, it can be seen that counter-gradient turbulent heat fluxes occur frequently deeper inside the tree layer in strongly stable conditions.

Similarity functions of wind and temperature
Figure 2g shows the profiles of the similarity function of wind (φ m ) in the three stability classes.On average, φ m has the largest value in the top layer of the canopy, being 3.2, 5.8, and 9.5 for the weakly, moderately, and strongly stable conditions, respectively.φ m decreases with depth in canopy and is usually negative in the middle layer in moderately and strongly stable conditions, due to the presence of the negative velocity gradient.In the lower layers, φ m values about 0.1 to 0.2.
In the top layer of the canopy, mean φ m increases with stability, different from the velocity gradient in Fig. 2b.This is because the local u * in the top layer is large (0.4 m s −1 ) in weakly stable conditions and decreases quickly with stability to 0.08 m s −1 in strong stratification.In the lower half canopy, the local u * is generally small and varies little with height, being 0.08 m s −1 in weakly stable conditions and about 0.05 m s −1 in moderately and strongly stable conditions.
The profiles of the similarity function of temperature (φ h ) in different stable conditions are presented in Fig. 2h.The profiles of φ h are different from those of φ m .In weakly stable conditions, φ h on average shows little sensitivity with height in the upper canopy but increases from about two to six with canopy depth increasing.However, in strongly stable conditions, the mean φ h is 43 at the top level, the largest in the profile, and decreases almost linearly with depth until to −17 in the layer just above the understory, with a value of 0.6 in the understory.In moderate stratification, larger variances of φ h are noted in the lowest two layers, a result of the very small values of covariance w θ v in those layers.

Relationship of the local stability parameter z/L(z) and Ri to the stability parameter h/L(h)
It has been shown in Sect.3.1 that thermal stratification exists within the forest canopy and that it is modulated by atmospheric stability conditions present above the canopy.Therefore, the atmospheric stability value determined above or at tree top may not be the best parameter to describe in-canopy stability.In Fig. 3  The two local stability parameters in the top layer (1.0-0.89h) have a good linear relationship with the stability parameter h/L(h) at the tree top on log-log scales, although local Richardson number Ri is more scattered with reference to h/L(h) than the mean local z/L(z) (Fig. 3a and d).The local Ri in the top layer linearly increases with h/L(h), without a horizontal asymptote predicted by Businger et al. (1971) (though it was originally for z/L < 1) as shown in Fig. 3d.This is different from that above canopy (or over the ground).Previous studies have shown that the gradient Richardson number Ri tends toward a constant (about 0.2) corresponding to large values of local stability parameter z/L (>1, up to 100) at 5.8 m above the Northern Spanish Plateau (Yagüe et al., 2006) and over snow and ice surface (Grachev et al., 2007).In the layer of 1.4-1.0h in our data set, local Ri also asymptotically tends to about 1.0 (higher than 0.2) when h/L(h) > 1 (data not shown).Ha et al. (2007) indicated, with CASES-99 data, that the gradient Richardson number Ri in weak wind, stable conditions is usually greater than predicted from z/L by Monin-Obukhov similarity theory (i.e., the asymptote of 0.2 in Ri for z/L > 1; Businger et al., 1971).This finding is consistent with our results in the top layer of the canopy where winds are usually much smaller than above the canopy.
At approximately mid-canopy (0.59-0.48 h) and in the lower layer (0.21-0.08 h) (Fig. 3b, c, e, and f), negative values of local stability parameters indicating local unstable conditions are often observed, especially when the atmosphere at tree top is very stable as indicated by h/L(h).When h/L(h) > 0.2, negative values of mean local z/L(z) occurs more frequently.Almost all cases are of negative values of mean local z/L(z) in the middle layer when h/L(h) > 1 and in the low layer when h/L(h) > 10 (Fig. 3b and c).Negative values of mean local z/L(z) in the middle and low layers reflect the positive covariance w θ or upward heat transport (Fig. 2f).The absolute values of the mean local z/L(z) tend to increase with stability h/L(h) in all layers on log-log scales, but are more scattered with depth inside the canopy (Fig. 3b  and c).
The local bulk Richardson number Ri in the middle and low layers versus h/L(h) is widely scattered with a large range of Ri values, especially when h/L(h) > 0.1, although there is a trend that Ri increases with h/L(h) when h/L(h) < 0.1 and decreases with h/L(h) when h/L(h) > 0.1 (Fig. 3e and f).Therefore, the two local stability parameters in the middle and low layers are independent of the stability parameter h/L(h) at tree top in stable atmospheric conditions.The latter cannot represent the local stability within the forest except in the upper layer of the crown.

Local flux-profile relationship for wind speed in the canopy layer
Local flux-profile relationships for wind speed (φ m ) within the canopy versus different atmospheric stability parameters are compared in Fig. 4. To be concise, only the data in the top (1.0-0.89h), middle (0.59-0.48 h) and low (0.21-0.08 h) layers are shown, corresponding to the layers of crown, trunk, and understory, respectively.With h/L(h) at tree top as the stability parameter, φ m in the top layer increases linearly with h/L(h) when h/L(h) < 1, which is shown as φ m first changes little (about 3) when h/L(h) < 0.1 and then increases fast on log-log scales (Fig. 4a).When h/L(h) > 1, φ m is so scattered that their relationship is less obvious.The variation of local φ m in the top layer with h/L(h) is quite similar to the relationships between φ m and z/L over 25-cm tall senescent grass (Howell and Sun, 1999), over snow and ice (Grachev et al., 2005(Grachev et al., , 2007) ) and over the ground (Yagüe et al., 2006), except that, when z/L is close to zero, φ m is close to 1 in their results comparing to about 3 in our results.In the middle and low layers of the forest canopy (Fig. 4b  and c), φ m values are far more scattered than those in the crown, with no obvious relationship with h/L(h).Negative φ m occurs especially in the middle layer when h/L(h) > 0.2-0.3,due to a negative wind speed gradient.Our results show that the relationship between the non-dimensional shear and h/L(h) agrees better with the Monin-Obukhov stability function in the crown than in the lower levels.This is consistent with the results in pine and aspen forests reported by Lee and Mahrt (2005), though there were no negative values of φ m included in their paper.
Using the mean local z/L(z) as the local stability parameter does not improve the flux-gradient relationship of φ m within the canopy (Fig. 4d-f), comparing to the relationship between φ m and h/L(h).In the top layer, the relationship between φ m and the mean local z/L(z) is very similar to that between φ m and h/L(h).In the middle and low layers, there is no obvious relationship between local φ m and mean local z/L(z) either, though there are two more combinations of positive/negative local φ m and negative mean local z/L(z).
The coefficient of power b is the slope with which local |φ m | linearly decreases with local |Ri| on log-log scales when |Ri| >1 (Fig. 4g-i).The slope is in the range of −0.41 to −0.47, averaged as −0.44 ± 0.023, for all layers (Table 2).This value is close to −0.5 that is the power of Ri when expressing local φ m as a function of local Ri by combining Eqs. ( 1) and ( 4).The coefficient b is not exactly equal to −0.5, which is because the other factors in the expression (e.g.u * and difference in temperature) are also related to Ri. Yagüe et al. (2006) showed the dependence of φ m on the gradient Richardson number Ri in stable conditions (up to Ri = 10) at different heights above the ground.The relationship between φ m and Ri they reported is similar to ours in www.biogeosciences.net/7/3625/2010/Biogeosciences, 7, 3625-3636, 2010 Fig. 4g for the top layer of the canopy, except that in their results φ m is about 1 to 2 when Ri < 0.02 and arrives the maximum values at Ri = 0.3.Mahrt (2008) also presented the dependence of φ m on Ri in stable conditions with Ri up to 10 3 above the ground.In his results, the relationship between φ m and Ri is similar to ours in the top layer (Fig. 4g) and above the forest (data not shown) when Ri < 1, but φ m shows little relation to Ri when Ri > 1, which is different from our results.Grachev et al. (2007) showed similar relationship between φ m and Ri but only for Ri < 1.In our results, φ m above the forest (1.4-1.0 h) is more similar to Yagüe et al. (2006) and also arrives the maximum values at Ri = 0.3, but tends to 2-3 when Ri < 0.02 (not shown).Within the canopy, local φ m keeps constant (about 3) for local Ri < 0.02 in the top layer and for local Ri < 0.1 in the middle layer.Previous studies also show constant φ m (about 1) for Ri < 0.02 over the ground (Yagüe et al., 2006;Grachev et al., 2007;Mahrt, 2008).

Local flux-profile relationship of temperature inside the canopy
Local flux-profile relationship for temperature (φ h ) is plotted against different atmospheric stability parameters throughout the canopy in Fig. 5.In the top layer of 1.0-0.89h, the relationship between φ h and the stability parameter h/L(h), local z/L(z), or local Ri is similar to each other.On log-log scale plots, φ h changes little with stability when h/L(h) < 0.1, mean z/L(z) < 0.1, or Ri < 0.02, but increases with stability when h/L(h) > 0.1, mean z/L(z) > 0.1, or Ri > 0.02, with more scattered φ h in strongly stable conditions (Fig. 5a, d,  and g).φ h decreases with h/L(h) and mean local z/L(z) when they are larger than 10 (Fig. 5a and d), but seems independent to Ri when Ri > 1 (Fig. 5g).These results in the top layer of the canopy are similar to the relationship of φ h to surface (or local) stability parameter z/L (or z/ ) over snow/ice surface by Grachev et al. (2005) and at the level of 5.8 m over the Northern Spanish Plateau by Yagüe et al. (2006), and to the bulk Richardson number Ri (Grachev et al., 2007;Yagüe et al., 2006).However, Grachev et al. (2007) also reported that φ h tends to be a constant for local z/L > 10 when data with temperature difference between the air and the snow surface less than 0.5 • C were omitted.
In the middle and deeper layers in the canopy, local |φ h | increases with h/L(h) when h/L(h) < 0.5, and decreases with h/L(h) when h/L(h) > 0.5 (Fig. 5b and c).With mean local z/L(z) as the atmospheric stability parameter, the observed magnitude of the mean local z/L(z) is usually larger than 0.1, and the local |φ h | decreases with the magnitude of the mean local z/L(z) (Fig. 5e and f).
With Ri as local stability parameter, when both local φ h and Ri are positive (blue circles in Fig. 5h and i) (i.e., both turbulent heat flux and temperature gradient are in downward direction), φ h still has a trend that φ h increases with Ri up to Ri = 1 and seems to be independent of Ri when Ri > 1 (Fig. 5h and i).This result is similar to that of Yagüe et al. (2006).They showed that, at different heights above the ground, φ h increases with Ri when Ri < 0.3 and then does not change significantly with Ri when Ri > 0.3.When local φ h and Ri have opposite signs, either positive Ri and negative φ h (red triangles in Fig. 5h and i) or negative Ri and positive φ h (green pluses in Fig. 5h and i), the points are much more scattered.These points correspond to either the combination of positive potential temperature gradient and upward heat flux or the combination of negative potential temperature gradient and downward heat flux, i.e., counter-gradient fluxes often seen in strongly stable conditions (Fig. 2e and  f).This result suggests that the presence of counter-gradient inside the canopy layer found in strongly stable conditions in the present study eliminates the otherwise valid local φ h and Ri relationships.

Self-correlation
We have seen from Sect.3.3 that, in the original data of φ m and Ri, their relationship is linear for |Ri| <1 and powerfunctional for |Ri| >1.One way to estimate the degree of self-correlation consists in first separating the original data into two groups with |Ri| <1 and with |Ri| >1 and then in separately examining the self-correlation test for each group.Correlating all the new points of one group together seems to eliminate the original relationship, because the randomly created Ri values may change from one range to another and follow different relationships between φ m and Ri.This will decrease the self-correlation and lead to possible misinterpretation.On the other hand, if truncating this new category of points for |Ri| < 1 (or |Ri| > 1), the analysis is no longer truly random (Klipp and Mahrt, 2004).Therefore, an alternate method is used to estimate the selfcorrelation as below.Each variable used to calculate φ m and Ri is picked at random from the pool of the entire original data set and new φ m and Ri are calculated.The selfcorrelation with |Ri| <1 (or with |Ri| >1) is estimated by calculating the linear correlation coefficient between new φ m and Ri for new |Ri| <1 (or by calculating the coefficient of determination of the nonlinear relation between new φ m and Ri for new |Ri| >1).The coefficient of determination for all Ri is also calculated with segment functions (i.e., linear function for new |Ri| <1 and nonlinear function for new |Ri| >1) to determine the self-correlation for all Ri.
The results of the self-correlation between φ m and Ri are listed in Table 3.It can be seen that, in the top layer (1.0-0.89h), the linear correlation coefficient (R) of the linear function with local |Ri| <1, the coefficient of determination (R 2 ) for the power function with local |Ri| >1, and R 2 for the segment functions for all local Ri of the original data set are 0.735, 0.744, and 0.494, respectively.They are much larger than the corresponding values of −0.139, 0.32, and 0.117 for random-created data set.This supports the idea that the original relationships between φ m and Ri in the top layer nighttime stable conditions have, at least partly, a physical meaning.
However, in the layers deeper inside the canopy, the original data set either shows |φ m | not changing with |Ri| for local |Ri| <1 in the middle layer (0.59-0.48 h), or only contains a few points with |Ri| <1 in the low layer (0.21-0.08 h).For local |Ri| >1, the coefficients of determination of the randomized data are 0.603 and 0.626 in the middle and low layers, very close to those (0.7 and 0.658) of original data set.Their difference is 0.098 in the middle layer and 0.032 in the low layer, much smaller than that in the top layer (0.424).This indicates large self-correlation for local |Ri| >1 in the two layers.A large self-correlation between local |φ m | and local |Ri| also exists in the two layers when the entire data set is taken into account (Table 3).
A similar procedure is applied to identify the degree of self-correlation between the local flux-profile relationship for temperature (φ h ) and the local bulk Richardson number Ri.Because φ h linearly increases with Ri for |Ri| <1 and tends to www.biogeosciences.net/7/3625/2010/Biogeosciences, 7, 3625-3636, 2010 For local |Ri| <1 just below tree top, the self-correlation coefficient is 0.28, a value comparable to the linear correlation coefficient (0.374) of the original data set with |Ri| <1, leading to small difference in R 2 (Table 4).The small correlation coefficient of the original data with |Ri| <1 is caused by the large φ h values of several points with Ri between 0.5 to 1.0 (Fig. 5g).When considering |Ri| <0.4, the correlation coefficient of original data is raised to 0.549 while the selfcorrelation coefficient keeps in the same level of 0.28.This makes the difference in R 2 go up to 0.221 (Table 4), indicating the physical meaning of the linear relation between φ h and Ri for small Ri in the top layer.It can also be seen from Table 4 that for |Ri| <0.4 in the middle layer, the linear correlation coefficient for local |φ h | and |Ri| is 0.374, larger than self-correlation coefficient of 0.181.This indicates that the relationship between φ h and Ri cannot be totally explained by self-correlation.Klipp and Mahrt (2004) showed that the gradient Richardson number Ri has less self-correlation with flux-profile relationships for wind speed φ m than the stability parameter z/L, where they used a different definition of φ m (≡ (kL/u * )(∂U /∂z)), in which the local Obukhov length "L" was used to replace "z" in Eq. (4).Mahrt (2007) examined the self-correlation between non-dimensional gradients (φ m and φ h ) and different stability parameters (z/L and Ri) in stable conditions with z/L < 1 or Ri < 0.25, above a grass surface.In that range of stability, their results support the presence of a physical relationship between φ m and Ri and between φ h and z/L, but not between the other two combinations of variables (φ m and z/L, φ h and Ri).We examined the self-correlation between the local non-dimensional gradi-ents (φ m and φ h ) and the local bulk Richadson number (Ri) within the canopy.The results show that, in the top layer of the canopy, the relationships (i.e., φ m linearly increasing with Ri for Ri < 1 and decreasing with Ri in a power function for Ri > 1 and φ h increasing with Ri for Ri < 0.4) are not totally caused by self-correlation and thus are physically meaningful.

Summary and conclusions
The relationships of non-dimensional wind and temperature gradients (φ m and φ h ) with different atmospheric stability parameters, i.e., h/L(h) at tree top, mean local z/L(z), and local bulk Richardson number Ri, within a tall forest canopy have been compared each other in stable atmospheric conditions.
The profiles of variables used in the calculation of φ m , φ h , and stability parameters suggest that, in strongly stable conditions, thermal stratification is generally present in the crown, trunks, and understory of the canopy.As expected, the accompanying velocity gradient is often small in the lower half canopy.In strongly stratified conditions, reverse velocity gradients often occur in the mid-canopy region with an accompanying upward transport of momentum observed in the lower canopy.During these periods, the presence of counter-gradient turbulent heat fluxes both in the middle and low canopy regions is also a frequent occurrence.
When the stability parameter h/L(h), the local |mean z/L(z)|, or the local |Ri| is less than 1, |φ m | increases linearly with them in the top layer (1-0.89h), but not in deeper layers in the canopy.When the magnitude of the three stability parameters is larger than 1, |φ m | has no clear relationship with h/L(h) and the local |mean z/L(z)|, but decreases with |Ri| in a power function and exhibits a general behavior throughout the forest.Therefore, local φ m has the best relationship with the local Richardson number Ri.
The local φ h in the top layer increases linearly with the two parameters h/L(h) and local |mean z/L(z)| when these are smaller than 1.When the local downward turbulent heat flux is consistent with the local temperature gradient in different layers of the canopy, the local φ h increases with the local Ri for Ri < 1 and does not change with Ri for Ri > 1.However, when counter-gradient heat transfer occurs inside a deeply stratified canopy, the relationship between the local φ h and Ri disappears.
The self-correlation analysis used here supports that the relationships between φ m and Ri (i.e., φ m linearly increasing with Ri for Ri < 1 and decreasing with Ri in a power function for Ri > 1) and between φ h and Ri (i.e., φ h increasing with with Ri for Ri < 0.4) at the top layer during nighttime are at least partly physically meaningful.

Fig. 1 .
Fig. 1.Vertical profile of leaf area index (LAI) density and its standard error of mean with sample number in parentheses, measured in June 2003 when the forest height (h) was 13.1 m.The total LAI was 2.85 ± 0.99.

Fig. 2 .
Fig. 2. Vertical profiles of (a) horizontal wind speed, (b) vertical gradient of horizontal wind speed, (c) covariance of horizontal and vertical velocities (u w ) normalized by friction velocity (u * ) at tree top, (d) difference of potential temperature (θ ) from that at tree top, (e) vertical gradient of θ , (f) covariance of w and θ normalized by u * and temperature scale (θ * ) at tree top, (g) similarity function of wind (φ m ), and (h) similarity function of temperature (φ h ), in weakly, moderately and strongly stable conditions with the atmospheric stability indicated with h/L(h) at tree top.

Fig. 3 .
Fig. 3. Comparison of local stability parameters, i.e., mean local z/L(z) and local Richardson number Ri in the top (1.0-0.89h), middle (0.59-0.48 h) and low (0.21-0.08 h) layers, with the atmospheric stability parameter h/L(h) at tree top.The dotted curve in (d) is the predicted Ri from h/L(h) (replacing z/L) with the formula given by Businger et al. (1971).

Fig. 4 .
Fig. 4. Comparison of the relationship of the local similarity function of wind φ m with different atmospheric stability parameters, i.e., h/L(h) at tree top, mean local z/L(z), and local bulk Richardson number Ri, in the top (1.0-0.89h), middle (0.59-0.48 h) and low (0.21-0.08 h) layers.The lines in (g-i) are the simulated relationships between |φ m | and |Ri| for |Ri| >1 on log-log scales.

Table 2 .
The simulated power function |φ m | = a|Ri| b with local |Ri| >1 in different layers within the canopy during nighttime, where φ m is the local flux-profile relationship of wind speed, Ri is the local bulk Richardson number, and a and b are coefficients.-0.08 162 0.683 −0.448 0.658 z/h, relative height of the layers, with h as the canopy height; N, the number of data points with local |Ri| >1; a and b, coefficients of the simulated function; and R 2 , the coefficient of determination.With the local bulk Richardson number Ri as the local stability parameter, when Ri < 1, local φ m increases linearly with Ri in the top layer (Fig. 4g), similar to with h/L(h).The increase of local φ m with Ri vanishes in the middle layer and local φ m does not change obviously with Ri when local |Ri| <1 (Fig. 4h).There are only a few cases with local |Ri| <1 observed deep in the canopy (Fig. 4i).When local |Ri| >1, local |φ m | decreases rapidly first and then slowly with the increase in local |Ri|, which is shown on log-log scales as linear decrease of local |φ m | with increase of local |Ri|, with a similar slope in all layers within the canopy (Fig. 4g-i).Local |φ m | also decreases with depth in the canopy for the same local |Ri|.The relationship between local |φ m | and local |Ri| with |Ri| >1 in different layers within the canopy can be described with a power function (Table 2), i.e., |φ m | = a|Ri| b ,for|Ri| > 1, (7) where a and b are coefficients of the power function.The coefficient a is the value of local |φ m | when local |Ri| = 1.It decreases with depth rapidly in the crown and then slowly in the lower layers, which can be described as an exponential function of layer height z/h:

Fig. 5 .
Fig. 5. Comparison of the relationship of the local similarity function of temperature φ h with different atmospheric stability parameters, h/L(h) at tree top, mean local z/L(z), and local Richardson number Ri, in the top (1.0-0.89h), middle (0.59-0.48 h) and low (0.21-0.08 h) layers.

Table 1 .
Description of the classification of stability based on the parameter h/L(h) at tree top, where h is the tree canopy height and L(h) is the Obukhov length at tree top.Due to removal of erroneous sonic anemometer data at 0.21 h in some days, the period numbers at the level 0.21 h or in the layers 0.48-0.21h and 0.21-0.08h are 15, 20, and 29 for weakly, moderately, and strongly stable conditions, respectively, less than those at upper levels, and thus the mean and standard deviation of h/L(h) are slightly different from those listed here.

Table 3 .
Self-correlation analysis between the local flux-profile relationship of wind speed, φ m , and local bulk Richardson number, Ri, in forest canopy, for a linear relation with |Ri| <1, a power function with |Ri| >1, and these segment functions for the entire data set./h, relative height of the layers, with h as the canopy height; N, the number of data points; R, correlation coefficient; R 2 , the coefficient of determination; subscript orig, for original data set; and subscript rand, for the average of 1000 trials of randomized data set. z

Table 4 .
Self-correlation analysis between the local flux-profile relationship of temperature, φ h , and local bulk Richardson number, Ri, in forest canopy, for linear relation with |Ri| <1 and |Ri| <0.4./h, relative height of the layers, with h as the canopy height; N, the number of data points; R, correlation coefficient; R 2 , the coefficient of determination; subscript orig, for original data set; and subscript rand, for the average of 1000 trials of randomized data set.beindependent of Ri for |Ri| >1 of the original data set near tree top (Sect.3.4),wefocus on the results of self-correlation test with small |Ri| listed in Table4. z