Is It Worthy to Apply Different Methods to Determine Latent Heat Fluxes? - A Study Case Over a Peach Orchard

This book represents an overview of the direct measurement techniques of evapotranspiration with related applications to the water use optimization in the agricultural practice and to the ecosystems study. Different measuring techniques at leaf level (porometry), plant-level (sap-flow, lysimetry) and agro-ecosystem level (Surface Renewal, Eddy Covariance, Multi layer BREB), are presented with detailed explanations and examples. For the optimization of the water use in agriculture, detailed measurements on transpiration demands of crops and different cultivars, as well as results of different irrigation schemes and techniques (i.e. subsurface drip) in semi-arid areas for open-field, greenhouse and potted grown plants are presented. Aspects on ET of crops in saline environments, effects of ET on groundwater quality in xeric environments as well as the application of ET to climatic classification are also depicted. The book provides an excellent overview for both, researchers and student,s who intend to address these issues.

energy balance equation as an indirect test of latent heat flux, LE, reliability (Brutsaert, 1988;Twine et al., 2000;Wilson et al., 2002;Foken et al, 2006;Castellví et al., 2008). The Penman-Monteith, PM, equation involves the aerodynamic resistance, but also the canopy resistance which makes PM equation difficult to apply for estimating LE in other than grazed surfaces not short of water (Allen et al., 1996). Sap flow measurements provide a method to measure transpiration, but despite it is not generally applicable, evaporation is not accounted which limits its application to irrigation system having little or no surface wetting. Soil moisture monitoring is widely used to estimate ET, but water movement between soil layers (even when the soil is unsaturated) is not easy to model, and soil moisture sensors measure water content within a small volume and the site measurement is often not representative of the entire field. Hence, due to the difficulty to obtain accurate ET (or LE) time series, even when affordability to direct measurement (i.e., using lysimeters and the EC method) is not a problem, simultaneous application of alternative methods or approaches may be especially useful when they do not share shortcomings. However, it is desirable to apply simple and affordable methods if they are proven reliable. An alternative method may help to perform quality control including gap filling. For LE estimation, the micrometeorological approach known as the residual method has been widely used (Brutsaert, 1988;Twine et al., 2000) Where, Rn is the net radiation, G is the soil heat flux and H is the sensible heat flux. In Eq. (1) the total rate of energy storage and other additional averaged energy sources (or sinks) including advective terms were neglected. However, Eq.
(1) appears to hold for most agricultural surfaces having adequate fetch , including tree orchards (Teixeira et al., 2008). When the available net surface energy, (Rn -G) is available, Eq. (1) combined with a method to estimate H allows for LE estimation. The objective of this study is to compare the hourly and daily series of LE determined over a peach orchard using a large weighing lysimeter, LE Lys , and from Eq. (1) estimating H using the EC method, LE EC , and the Surface renewal, SR, method, LE SR Castellví, 2004;Paw U et al., 2005). The SR approach was selected as an alternative method because the instruments required are affordable, robust, easy to maintain, transportable, sense much larger area than a lysimeter, and it is expected reliable under windy conditions. For this study, the SR model used combines SR analysis, mixing-length theory for momentum (Harman and Finnigan, 2007) and mixing-layer analogy (Raupach et al., 1996;Graefe, 2004). Therefore, it operates close to the canopy top which avoids installation of tall mast as required for other micrometeorological methods. To preserve stationarity, hourly and daily LE estimates (Eq. 1) were determined integrating the half-hourly values of the measured available net surface energy and sensible heat flux. The LE estimates were compared with measurements from a weighing lysimeter located nearby.

Theory
SR analysis [pioneered by Paw U et al. (1995)] combines SR theory [pioneered by Higbie (1935)] and the analysis of the scalar time trace (a subject of research during decades that still is of major interest). The SR theory (originally developed to investigate interfacial heat transfer between a liquid and a gas) assumes that heat transport occurs when a fluid parcel travelling at a given height in the bulk of the flow (above the interface) descends and takes contact with the heated surface for a period during which the parcel is heated. There is an www.intechopen.com unsteady diffusion transport during the contact until the parcel is, by continuity, replaced or renewed by other parcel coming from aloft. The detection of sequential sweeps and ejections of parcels from the interface is crucial because the heating take place in between.  applied the SR concept over land surfaces to investigate scalar transfer between sources located at the surface and the atmospheric surface layer. That is, consider an air parcel travelling at a given height above the canopy. Due to shear stress, SR analysis assumes that at some instant the parcel moves down into the canopy and remains in contact with the canopy elements for a period after it is ejected upwards and replaced by another parcel sweeping in from aloft. During the connect period the parcel has been enriched (depleted) of scalar due to the exchange between the air and sources (sinks). Each sweep and ejection is an injection of scalar from the sources into the bulk of the atmosphere which is characterized by a given amount of scalar released during the renewal period. The total flux of a scalar is, therefore, driven by the continuous renewal (i.e., replacement of air parcels) across the canopy top averaged for a given period (typically, half-hour). For sensible heat flux, when high-frequency air temperature measurements are taken at a given height z, the renewal process can be visually inferred in the time trace as a rather regular low-frequency ramp-like (asymmetric triangle shape) pattern. Paw U et al. (1995) presented a diagram of the SR process ( Fig. 1a) and abstracted an ideal scheme for a ramp-like event in the air temperature trace (scheme 1 in Fig. 1a). Chen et al. (1997a) presented a slightly different ramp model (scheme 2 in Fig. 1a) that neglects the quiescent period but includes a micro-front period instead of an instantaneous ejection. Whatever the scheme, the amplitude, A, and period, , are the ramp dimensions that characterize an injection of sensible heat flux. The eddy responsible of such injections of scalar (i.e., the coherent structure) explains most of the total flux determined at the measurement height using the EC method (Hongyan et al., 2004). Analysis of time series, therefore, consists on identification of coherent structures to extract their ramp dimensions as shown in Fig. 1b. A method based on structure functions for determining the mean ramp dimensions sequentially observed in a trace (typically half-hour) is described in Appendix A. Based on the scalar conservation equation for a planar homogeneous turbulent flow, assuming that the air parcel is uniformly heated by sources (sinks) below the measurement height, z, with no heat lost from the parcel top while it remains in contact with the canopy, the mean sensible heat flux density can be estimated as  () where and C p are the density and specific heat at constant pressure and, z, is the top of the air parcel which represents a volume per unit area. The parameter  is included to correct for the assumptions invoked. When measurements are taken above the canopy, based on the one dimensional diffusion equation, the following relationship to estimate over the averaging period was derived (Castellví, 2004)   The time course of the air temperature within the air parcel for the different positions is idealized in two ramp models. Scheme 1 assumes a quiescent period and a sharp instantaneous drop in temperature. Scheme 2 neglects the quiescent period and assumes a finite micro-front. L r , L q and L f denote the warming, quiescent and micro-front periods, respectively. A is the ramp amplitude and  is the total ramp duration. B) Air temperature measured at 10 Hz versus time observed under unstable conditions over a peach tree during 30 s. Few ramps fitted (by eye) are shown for the first 5 s.
Equations (2) and (3) provide a method exempt of calibration for estimating sensible heat flux (Snyder et al., 1996;Castellví et al., 2008). In Eq. (3), d is the zero-plane displacement, z* is the roughness layer depth (height from the ground to the base of the inertial sub-layer), h is the canopy height, u * is the friction velocity, k = 0.4 is the Von Kármán constant,  h ( is the stability function for heat transfer described below in Eq. (5), and  is an stability where g is the acceleration due to gravity and T v is the virtual temperature, which can be replaced with T in dry environments. Likely, the most widely used formulation for  h () is (Högström, 1988;Foken, 2006) It is of interest to mention that other expressions covering a wider range of stability conditions are available (Kader and Perepelkin, 1989). www.intechopen.com

Procedure to estimate the sensible heat flux
To solve Eq.
(2), other equations providing estimates for z*, d and u * are required. Therefore, an iterative procedure to solve the sensible heat flux must be implemented. Certainly, to achieve the full potential of SR analysis, the input required should be constrained to measurements taken at a single height. A model described in Castellví and Snyder (2009) which is based on a mixing-length theory for momentum (Harman and Finnigan, 2007) combined with mixing-layer analogy (Raupach et al., 1996) was used for half-hourly z* estimates as a function of the canopy architecture (canopy height, h, leaf area index, LAI, and the crown thickness relative to the canopy height, f), drag coefficient at leave scale, c d , and the turbulent intensity near the canopy top, I u (= u h u  ; where u is the turbulent standard deviation of the horizontal wind speed, u, and u h is the mean u at the canopy top). It is described in Appendix B. Appendix C shows the procedure used to solve H. Accordingly, z* estimation taking measurements taken at one height may be obtained using the following expression

The field experiment and main climate features
The experiment was carried out at the University of California Kearney Research and Extension Center in Parlier (CA) from August 2 nd to October 16 th 2007. The peach orchard (Crimson Lady) architecture consisted on mature trees, 3.95 m tall, with a dense crown from 0.75 to the canopy top, 4 m distance between trunks in a row and 4.5 m between rows. The leaf area index was LAI≈3.0 and the ratio of the depth of foliage (m) to the canopy height was f≈0.95. The orchard was mainly surrounded by grapevines and bare soil. The lysimeter consisted of an underground chamber that houses a balance-beam weighing system (Fred Lourence, Precision Lysimeters, Red Bluff, CA) which is (2m x 4m x 2m deep) containing two trees with similar spatial separation of those in the orchard. Weight changes were logged hourly and the error in latent heat flux was 0.05 mm h -1 (Scott et al., 2005). The lysimeter fetch in the prevailing wind direction was 240 m. As a rule of thumb, fetch requirements are estimated according to the ratio 1:100 (i.e., the adjusted surface sublayer grows 1 m for each 100 m distance to the leading edge in the wind direction) and the basis of the inertial sublayer is located of about two-three times the canopy height (Brutsaert, 1988;Kaimal and Finnigan, 1994). The latter prevents the use of traditional micrometeorological methods to estimate turbulent surface fluxes close to the lysimeter, especially those that require measurement of gradients. Latent heat flux, Eq. (1), was evaluated half-hourly nearby the lysimeter. To prevent potential differences due to influence of local advection, instrumentation above the canopy was deployed 40 m apart in the cross mean streamwise direction having same fetch. Net radiation was measured using a net radiometer (REBS, Inc Q7.1, Bellevue, Wa) placed at 6.0 m. The soil heat flux was measured as described in Fuchs and Tanner (1967) and it was obtained as the average of three measurements to account for spatial variability. The three wind speed components and sonic temperature were measured using a sonic anemometer (81000RE, RM Young, USA). It was deployed at 5.5 m to ensure measurements within the adjusted surface layer and to avoid the region of the flow with maximum absorption of momentum. The raw sonic data was recorded at 10 Hz using a data logger (CR1000, Campbell Sci., USA). The protocol used as a reference for comparison described in Mauder et al. (2007) was applied using the TK2 package [free distributed by the University of Bayreuth (Mauder and Foken, 2004)], to determine half-hourly means, variances and covariances. Half-hourly fluxes were used to evaluate the hourly fluxes in Eq.
(1). Because the latent heat flux estimates were determined as a residual of the surface energy balance (Eq. 1), the work of expansion of moist air parcels under constant pressure were included in the sensible heat flux. Therefore, the covariance of the vertical wind speed with the sonic temperature, which is close to the virtual temperature, was used because it mostly accounts for the work of expansion of the moist air parcel (Paw U et al., 2000). In practice, when air temperature is measured the work of expansion can be estimated as 0.076 times the latent heat flux (Paw U et al., 2000). Climate features. Clear skies and no rainfall were observed during the experiment which is the typical climate features in the San Joaquin valley. The following half-hourly means and standard deviations, respectively, were observed; 0.8 m s -1 and 0.5 m s -1 for the horizontal wind speed, 22.4 C and 7.1 C for sonic air temperature, 0.21 m s -1 and 0.15 m s -1 for the friction velocity, and 24 W m -2 and 70 W m -2 for H determined using the EC system; 106 W m -2 and 138 W m -2 for LE measured using the lysimeter; and 129 W m -2 and 196 W m -2 for (Rn-G). The number of half-hourly samples collected for unstable and stable cases was 1692 and 1904, respectively. The larger dataset for stable cases was due to the formation of a capping inversion during the afternoon due to regional advection of sensible heat flux in the San Joaquin Valley. Typically, the surface boundary layer becomes near neutral and stable from about two-three hours after noon until sunrise.

Results and discussion
To prevent stationarity, sensible and latent heat fluxes using the EC method and SR analysis were estimated half-hourly. The hourly and daily flux estimates were determined by integrating the half-hourly values and were denoted using the corresponding subscripts EC and SR. According to Eq. (1), LE EC and LE SR were determined as; Rn -G -H EC and Rn -G -H SR , respectively. To be consistent with the lysimeter (the reference), the fluxes were expressed in mm. Table 1 shows the linear regression analysis (slope, a, intercept, b, determination coefficient, R 2 ) and the root mean square error, Rmse, to compare the hourly H SR versus H EC for unstable and stable cases and for all the data. It is of interest to mention that, despite in micrometeorology the EC method is considered a reference it operates best when deployed well above the canopy top (i.e., in the inertial sublayer) because large eddies are easier to sample. Thus, likely the EC method slightly underestimated the actual sensible heat flux. Regardless, realistic H estimates must be well correlated with H EC because the EC system directly measures the turbulence. Table 1 shows that, regardless of the stability cases, the slopes were rather close to one, the correlations were high, and the intercepts and the Rmse were small. The performance between H SR and H EC is shown in Fig 2A and, subsequently, LE SR and LE EC were similar (not shown) regardless of the integration period (hourly and daily). Figure 2B compares the hourly LE EC versus LE Lys . Visually, differences observed between LE EC and LE Lys are in the order of H ( Fig. 2A). One may presume that integration of errors due to missing energy terms in the surface energy balance, errors in measuring the available net surface energy, deployment of the EC system to close to the canopy, natural spatial variability in fluxes over a heterogeneous surface enhanced from trees within a lysimeter and differences in footprints may explain most of the scatter observed in Fig. 2B. However, Fig. 3   not matched. Partially, it could be expected as a result of the scatter shown in Fig. 2B (i.e., errors did not balanced at daily scale). In Fig. 3 it is shown that during the first month, LE EC and LE SR were higher than LE Lys . During the period from day of the year 273 to 290, LE EC www.intechopen.com and LE SR fluctuated at a rate smaller than LE Lys . After looking at the lysimeter management reports, the different time course shown in Fig. 3 obeyed to irrigation problems. Blockage of the drip irrigation system led to insufficient water inside the lysimeter during the first month, so the ET from the lysimeter did not increase with evaporative demand like the remainder of the orchard, which had adequate soil moisture. As a consequence, after day of year 273 the lysimeter was flood irrigated to return the soil to well-watered conditions. Thus, the deficit irrigation due to the system problem explains the lower LE Lys during the first month and flooding of the lysimeter explains the higher LE Lys after day 273. For the campaign, the total ET estimated by LE EC was 4.8% lower than LE SR , and both were close to the total LE Lys (i.e., within 2%). Partly, this issue can be explained because the extra amount of water applied to restore the lysimeter was evaluated by the difference in the irrigation observed in and out of the lysimeter.  (Castellví and Martinez-Cob, 2005;Castellví et al., 2006a). Therefore, appears that at sites were light winds are often observed application of SR analysis requires estimation of z*. Though not directly comparable because measurements were taken in the inertial sublayer and fetch was large, a study carried out over short drip irrigated grass (0.10 -0.15 m tall) under similar weather features, comparison of hourly LE EC and LE SR versus LE Lys showed excellent agreement (Castellví and Snyder, 2010). The latter was as a consequence that, despite over irrigated surfaces (Rn -G) explains much of LE, adding in Eq. (1) the H contribution the correlation of the hourly LE EC and LE SR versus LE Lys was significantly improved.

Summary and concluding remarks
An experiment was carried out in a peach orchard to study the reliability of the latent heat flux estimates using the EC method and SR analysis for estimating sensible heat flux in conjunction with the surface energy balance equation. Therefore, it is assumed that the net surface energy (Rn -G) is available. The fetch was limited, the site is influenced by regional advection of sensible heat flux, and light winds occur during most of the day. As a consequence, measurements were taken close to the canopy and the roughness sub-layer depth was estimated on a half-hourly basis to better estimate the  parameter required in SR analysis. A model based on a mixing-length theory for momentum combined with mixinglayer analogy allowed z* estimates from some characteristic canopy parameters and the turbulent intensity measured near the canopy top. To test the reliability of the z* estimates, detailed profiles of the wind speed and air temperature are needed. They were unavailable, regardless, determination of z* appears difficult due to limited fetch. It is shown the potential of SR analysis to estimate the sensible heat flux. It was found that SR analysis was highly correlated with the EC method. Latent heat fluxes estimated as, LE = Rn -G -H, were compared with lysimeter data, and the comparison clearly showed the impact of irrigation management on the lysimeter ET. Further research is required covering long campaigns, however, this study suggests that SR analysis for estimating LE using the residual method can be taken in consideration as methods that may help to perform quality control of LE time series.

Acknowledgements
The author thanks to R.L Snyder (Univ. of Calif., Davis, CA), N.Rambo and S. Ewert (Dept. of Water Res., San Joaquin District, Fresno, CA), S. Johnson (Univ. of Calif., Kearney Agric. Res. and Ext. Center, Parlier, CA), J. Ayars ( USDA ARS, Fresno, CA) and A. Russo (Univ. of Catania, IT) who carried out the field experiment, and to Asun, Carla and Tania for their help in using various facilities at Lleida (Spain). This work was supported by a grant (#4600004549) from the California Department of Water Resources (Sacramento, CA), Ministerio de Ciencia y Innovación and TRANSCLA project (CGL2009-12797-C03-01) of Spain.

Appendix A: Determination of the mean ramp dimensions using structure functions
It is of interest to mention that the technique based on structure functions of different order to determine the ramp dimensions is objective (i.e., there is not need to implement filters as for other methods), all the measured data is used, and assumes ramps in a sequence. Based on ramp model shown in scheme 1, Van Atta (1977) pioneered this technique that assumes a statistical independency between the coherent structure and the smallest eddies. The ramp model shown in scheme 2 is described. It accounts for a micro-front period. Therefore it is more realistic than the model shown in scheme 1 (Chen et al., 1997a). A structure function of order n is defined as where m is the number of data points in the 30-minute interval measured at frequency (f) in Hz, n is the power of the function, j is a sample lag between data points corresponding to a time lag (r = j/f), and T i is the ith temperature sample. Ramp dimensions A and t f can be determined by fitting the 2 nd , 3 th and 5 th order structure functions to the following equations for different time lags If t f is neglected in the above expressions, the resulting expression for 0 < r ≤  can be used to determine the ramp dimensions A and  according to scheme 1 (Van Atta, 1977). The mean ramp amplitude is determined by solving the following expression for the real roots:  Chen et al. (1997a), the shortest time lag to be used for linearization, r 1G , is that produces the first global maximum of S 3 (r)/r because A.1 does not hold for r < r 1G . Thus, the initial time lag, r ini , used to linearize A.1, was r ini =r 1G . To estimate the maximum time lag, r end t o b e u s e d f o r linearization so that r << Lr, the second global maximum of S 3 (r)/r was determined. Based on the third order structure function for , t f < r ≤ (t f ), the second global maximum occurs at a time lag, r 2G , giving r 2G ≈ ¾ . According to Qiu et al. (1995), Lq ≈ 0.25  in scheme 1, and therefore, r 2G ≈ Lr. The last r used to linearize A.1 was determined as 2% of r 2G to insure r end << Lr.

Appendix B: Estimating the roughness layer depth
Turbulence in the roughness layer is characterized by the presence of distinct coherent structures generated near the canopy top (Finnigan and Shaw 2000;Shaw et al. 2006). Therefore, similarity does not apply within the roughness layer (Kaimal and Finnigan, 1994). As a result many of the characteristic properties of the roughness sub-layer differ from those of the inertial sub-layer and resemble more those of a plane mixing layer [i.e., the mixing-layer analogy reported by Raupach et al. (1996)]. In the following, to simplify, it is assumed that the canopy is homogeneous and dense at an extent that is capable to absorb all the momentum transferred into the canopy by coherent structures. Therefore, because the location of the ground is irrelevant to the dynamics involved, the natural choice of z-axis origin is the location where the physical processes change. That is, the z-axis origin is set at the canopy top which is at height h. Accordingly, (z*h) defines the sub-layer that extends from the canopy top to the bottom of the inertial sub-layer. Within the canopy, the horizontal mean wind speed, u, and shear reach a maximum at the canopy top. They are attenuated within the canopy at a rate determined by the drag exerted by the canopy elements. Attenuation is mainly important in the upper part of the canopy and the vertical profile of the mean wind speed follows an exponential decay which can be expressed as (Kaimal and Finnigan, 1994)  where LAI is the leaf area index and f the crown thickness relative to the canopy height. When measurements are taken close to the top of the canopy top, the relationship u*  0.5 u holds (Kaimal and Finnigan, 1994), where  u is the mean u turbulent standard deviation.
Thus, we arrive to the approximation:

Appendix C: Determination of the sensible heat flux
Estimation of z* (B.8) presumes that measurements are taken at z=h. However, in practice instrumentation is deployed slightly higher to avoid potential shortcomings and damage. In general, depending on the stability conditions, the measurement height may fall within the roughness or in the inertial sublayer. Therefore, the correct expression to estimate  must be selected, and the friction velocity may be estimated as follows (Kaimal and Finnigan, 1994)  where, u r is the wind speed at reference height z r , and  m () is the integrated Businger-Dyer relationship for momentum (Paulson, 1970) 2 2 ln(0.5(1 )) ln(0.5(1 )) 2 arctan( ) 0.5 0 () 4.7 0 where, x is, x=(1-16) 1/4 . It is presumed that the wind profile is unavailable. Therefore, for measurements taken in the inertial sub-layer over dense canopies the zero-plane displacement and the aerodynamic surface roughness length, z o , can be estimated as a portion of the canopy height, h, d=0.67 h and z o =0.12 h (Brutsaert, 1988;Wieringa, 1993). The stability parameter and the ramp amplitude for temperature have different sign (Fig.1). Therefore, after determining the ramp dimensions (Appendix A), the appropriate expressions for  h () and  m () are known. Estimation of the mean drag coefficient at the leaf scale is a compromise. It depends on the shape and orientation of the leaves and, in principle, it may depend on the velocity field within the canopy through the Reynolds number though such relationship is not still clear (Brutsaert, 1988). Numerical adjustment to produce the best agreement between models of transfer of momentum within the canopy and observations is used to determine c d values (Ionue, 1981;Pingtong and Takahashi, 2000;Mohan and Tiwari, 2004). Because c d values are mostly indirectly determined under neutral conditions, if possible, it is best to carry out a short campaign for adjustment. Regardless of the stability conditions, the c d was set to, c d =0.2, which is a typical value (Graefe, 2004). It is of interest to mention that because the zero-plane displacement, roughness length and the friction velocity were estimated, the c d value obtained accounts for such uncertainties. A recent study carried out over orange trees covering a year of measurements shows that a short dataset of about two weeks is enough to adjust B.8 by trial and error to fit H SR and H EC .
Therefore, the SR model appears robust. Next, the sensible heat flux, , and u * are solved simultaneously by iteration. To start the iteration procedure, neutral conditions are assumed for the actual atmospheric surface layer. This gives an approximation of the actual friction velocity (case z > z*), parameter  and H. The latter values allow for the first approximation of and the process is iterated until convergence is achieved.