Calibration of METRIC Model to Estimate Energy Balance over a Drip-Irrigated Apple Orchard

A field experiment was carried out to calibrate and evaluate the METRIC (Mapping EvapoTranspiration at high Resolution Internalized with Calibration) model for estimating the spatial and temporal variability of instantaneous net radiation (Rni), soil heat flux (Gi), sensible heat flux (Hi), and latent heat flux (LEi) over a drip-irrigated apple (Malus domestica cv. Pink Lady) orchard located in the Pelarco valley, Maule Region, Chile (35◦25′20′′LS; 71◦23′57′′LW; 189 m.a.s.l.). The study was conducted in a plot of 5.5 hectares using 20 satellite images (Landsat 7 ETM+) acquired on clear sky days during three growing seasons (2012/2013, 2013/2014 and 2014/2015). Specific sub-models to estimate Gi, leaf area index (LAI) and aerodynamic roughness length for momentum transfer (Zom) were calibrated for the apple orchard as an improvement to the standard METRIC model. The performance of the METRIC model was evaluated at the time of satellite overpass using measurements of Hi and LEi obtained from an eddy correlation system. In addition, estimated values of Rni, Gi and LAI were compared with ground-truth measurements from a four-way net radiometer, soil heat flux plates and plant canopy analyzer, respectively. Validation indicated that LAI, Zom and Gi were estimated using the calibrated functions with errors of +2%, +6% and +3% while those were computed using the standard functions with error of +59%, +83%, and +12%, respectively. In addition, METRIC using the calibrated functions estimated Hi and LEi with error of +5% and +16%, while using the original functions estimated Hi and LEi with error of +29% and +26%, respectively.


Introduction
Water scarcity has become a main limiting factor for agriculture, especially in arid and semi-arid conditions.This issue has become particularly relevant in countries like Chile where many commercial apple orchards have been transformed from traditional to intensive (1500 tree•ha −1 ) or super-intensive (3500 tree•ha −1 ) planting densities that require the optimization of water application and reduction of the cost of electrical energy required for pumping water into orchards.For this reason, decision-making tools for irrigation management are needed to optimize water productivity and maintain sufficient levels of crop productivity and quality.A correct quantification of apple water consumption or actual evapotranspiration (ETa) is a key factor in improving water productivity both at orchard and regional scales.
Traditionally, ETa is calculated by multiplying the reference evapotranspiration (ETr) by a crop coefficient (Kc) which depends on crop type and the crop growth stage [1,2].This procedure is based on the Penman-Monteith (PM) equation, which computes daily ETr values over a reference surface (alfalfa or grass) using weather-based data from an automatic weather station (AWS) [3,4].However, the Kc values taken from literature are often not adapted to varying orchard conditions such as soil type, cultivar and climate characteristics [5,6].In addition, Kc values for complex, heterogeneous canopies such as fruit orchards depend highly on canopy architecture, plant densities, and standard agronomical practices that determine the partitioning of ETa into transpiration and soil evaporation [7].This issue can produce substantial error in the estimation of orchard water requirements and reduce the water use efficiency, yield and fruit quality [8].
In recent decades, remote sensing energy balance (RSEB) approaches have become a valuable technology to quantify the spatio-temporal variability of water requirements for a wide range of vegetation types [9][10][11].To map ETa from satellite imagery, RSEB algorithms use the remotely sensed surface reflectance and radiometric surface temperature (Ts) to transform satellite radiances into land surface characteristics such as emissivity, surface albedo, vegetation indices, and surface aerodynamic temperature [12][13][14].Thereby, satellite-based data in combination with ground-based meteorological measurements are used to compute the daily ETa as a residual from the surface energy balance equation at the time of the satellite overpass [15]: where λ is the latent heat of vaporization (J•kg −1 ); ETa i is instantaneous evapotranspiration (mm•h −1 ); LE i is instantaneous latent heat flux (W•m −2 ); Rn i is instantaneous net radiation (W•m −2 ); G i is instantaneous soil heat flux (W•m −2 ); and H i is instantaneous sensible heat flux (W•m −2 ).Subscripts "i" denote measurements at the time of the satellite overpass (approximately 11:30 h local time).RSEB models estimate Rn i , G i and H i using surface reflectance and surface temperature from satellite platforms such as MODIS, Landsat (7 ETM+ and 8 OLI) and ASTER [12,[16][17][18][19]. RSEB models have been tested extensively in different parts of the world with positive results for several fully covered crops but little research exists on the application of these algorithms to estimate ETa over heterogeneous or complex canopies such as apple orchards.In this regard, studies have indicated that the METRIC (Mapping EvapoTranspiration at high Resolution with Internalized Calibration) model is able to estimate ETa with errors ranging between 5% and 20% in fully covered crops such as beans, wheat, corn, cotton and beets [6,19,20].Several researchers have indicated that METRIC can be used to evaluate the spatial variability of ETa and Kc of vineyards and orchards with promising results [16,18,[21][22][23].
METRIC is a one-source model that uses an internal self-calibration process called CIMEC (Calibration using Inverse Modeling at Extreme Conditions) [3,6,12,16], which eliminates the need for an atmospheric correction of surface temperature and albedo measurements using the radiative transfer model [16,19].In addition, METRIC utilizes ground-based ETr data acquired in the satellite scene to inversely and internally calibrate the surface energy balance which facilitates determining the extreme conditions of heat exchange over agricultural areas (the "hot" and "cold" pixels or anchor pixels) in order to calibrate the sensible heat flux at the pixel-level [16,24,25].
The main disadvantage of the energy balance approach employed in METRIC is that the estimation of LE i for each pixel is only as accurate as the estimates for Rn i , G i and H i [6].In addition, sensitivity analyses applied to METRIC indicate that variations on intermediate parameters such as albedo, leaf area index (LAI), G i and momentum roughness length (Z om ) can indirectly affect the accuracy of LE i and subsequently ETa estimation [25,26].Therefore, the evaluation of energy balance over heterogeneous canopies such as orchards is needed to provide insight into sources of error in the estimation of orchard ETa.
In general, sparse woody canopies such as high-density apple orchards have spacings between rows that vary with canopy geometry and planting density and that affect the partitioning of Rn i into G i , H i and LE i [27].Therefore, a satellite-pixel integrates a mixture of different surface temperatures, e.g. the shaded and illuminated portion of the canopy in combination with the shaded and sunlit surface of bare soil at the moment of satellite overpass which can lead to bias in the estimation of H i [23] due to variation in sun angle and satellite viewing angle.To minimize these effects, a correction is required in the sub-models that compute G i , LAI and Z om [18].In addition, Santos et al. (2012) [23] reported that the estimations of H i may be biased by uncertainty in the modeling of Z om which is a measure of resistance and friction between the air layer that interacts with the surface of crops [17].In the case of fruit orchards, Z om increases in proportion to the crop density until a density threshold of approximately LAI = 3 m 2 •m -2 is achieved when the trees have a uniform distribution of the leaves in the canopy [3].
Therefore, the application of METRIC to estimate orchard water requirements requires specific adjustment of sub-models that estimate LAI, G i and Z om .Thus, the main goals of this research were to calibrate the sub-models to compute LAI, G i and Z om and to validate the METRIC model to compute the instantaneous values of Rn i , G i , H i and LE i over a drip-irrigated apple orchard under Mediterranean semi-arid climatic conditions.

Original METRIC Model
METRIC uses TIR-Multispectral satellite images (e.g., from Landsat 7 ETM+ and 8 OLI) and ground-based meteorological data to estimate pixel-based daily ETa at 30 × 30 m spatial resolution via a surface energy balance equation at the time of the satellite overpass (Equation ( 1)).Rn i (W•m −2 ) is calculated pixel by pixel as [17]: where α is surface albedo (dimensionless), R s↓ is the incoming short-wave radiation (W•m −2 ), R L↓ is the incoming long-wave radiation (W•m −2 ), and R L↑ is the outgoing long-wave radiation (W m −2 ).ε o is broadband surface thermal emissivity (dimensionless).The subscript " M " indicates values computed by METRIC.
In the METRIC algorithm [6], values for G i are estimated as follows [28]: where Ts is the radiometric surface temperature ( • C) estimated from the TIR band and NDVI is the normalized difference vegetation index.For homogenous canopies, values of c 1 and c 2 are 0.0038 and 0.0074, respectively.For short annual crops, LAI (m 2 •m −2 ) is estimated as [17]: where SAVI is the soil adjusted vegetation Index (dimensionless) calculated using the top of atmospheric reflectance.
Values for H i (W•m −2 ) are obtained for each pixel as [16]: where ρ air is air density (kg•m −3 ); C p is specific heat of air at constant pressure (1004 J•kg −1 •K −1 ); and ∆Ts is the near-surface air temperature gradient (∆Ta = Ta z1 − Ta z2 ) above each pixel, where Ta z1 and Ta z2 are near surface air temperature ( • K) at heights z 1 and z 2 above the elevation of d + Z om , where d is zero plane displacement height (all in m).r ah is the aerodynamic resistance to heat transport (s•m −1 ).Details on estimating r ah and H i_M are described by [16] and include correction for effects of buoyancy.
To estimate the ∆Ts, two anchor pixels are selected in the study area corresponding to cold and hot pixels in order to apply the CIMEC process [6,16].The cold pixel is defined and located over a well-watered and non-stressed agricultural area, which represents the maximum value of ETa, while the hot pixel is located in a bare soil agricultural area where ETa is nearly zero [29].The internal calibration reduces impact of any biases in estimation of aerodynamic stability correction and surface roughness [16].∆Ts is linearly calculated at all pixels as: where a and b are empirically determined constants for a given satellite image and Ts is the surface temperature ( • K) adjusted to a common elevation datum for each image/pixel using a digital elevation model and customized de-lapsing rate.METRIC estimates Z om of annual agricultural crops assuming that the crop height varies proportionally with the LAI, as [17]: Z om = 0.018•LAI (8)

Study Area
The study was conducted in 5.5 ha of drip-irrigated apple (Malus domestica cv.Pink Lady) orchard located in the Pelarco valley, Maule Region, Chile (35 • 25 2 LS; 71 • 23 57 LW; 189 m above sea level) (Figure 1) during three growing seasons (2012/2013, 2013/2014 and 2014/2015).The experimental plot was surrounded by other apple fields with similar management and conditions.The climate of the study area corresponds to a typical Mediterranean with a daily average temperature of 14.5 • C between September and March (local spring to summer) and average annual rainfall of 700 mm, concentrated mainly throughout the winter months.The summer period is usually hot (30 • C max.) and dry (3.5% of annual rainfall) with a very high atmospheric demand for water vapor, while the winter is cold (−2 • C min) and wet (74.5% of annual rainfall).The soil at the experimental field is classified as "Pelarco series", which corresponds to a loamy clay texture on an alluvial terrace position having flat topography (34% clay, 30% silt and 36% sand).For the effective rooting depth (0-90 cm), the average volumetric soil water content at field capacity (θ FC ) and wilting point (θ WP ) were 0.22 m 3 •m −3 (132 mm) and 0.35 m 3 •m −3 (210 mm), respectively.
The apple trees were planted in 2008 with a density of 1667 trees•ha −1 (4.0 × 1.5 m plant spacing) and irrigated daily using drippers (4 L•h −1 ) spaced at intervals of 0.75 m (Figure 1a).The experimental plot was maintained under non-water stress conditions during the whole season.The trees were trained on a vertical axis system with an average tree height of 4 m and canopy width of 1.55 m.Typical tree trunk diameters ranged between 65 and 79 mm during the three study seasons.

Plant Measurements
To evaluate the quality of irrigation management, the midday stem water potential (Ψ md ) was evaluated weekly with a pressure chamber (PMS 600, PMS Instrument Company, Corvallis, OR, USA) during the three growing seasons.Fully expanded leaves (two per tree, one tree per replicate) were removed from 18 trees, wrapped in aluminum foil, and encased in a plastic bag at least 2 h before measurement [30].Measurements were done at midday (12:00-14:00 h) from sunlit leaves situated on both sides of the apple rows.Volumetric soil water content (θ i ) at the rooting depth (0-0.6 m) was monitored weekly at 10 sampling points distributed inside the experimental plot using a portable time domain reflectometry (TDR) unit (TRASE, Soil Moisture Corp., Santa Barbara, CA, USA).LAI was measured biweekly in 20 trees during the study period using a plant canopy analyzer (LAI-2000, LI-COR, Lincoln, NE, USA) with diffuse radiation at sunset.Finally, the plant canopy analyzer was calibrated defoliating 4 trees, where the leaves of each tree were photographed, and then, the total leaf area per tree was measured using a digital image analysis [7,31].

Weather Data for ETr-PM Estimation
An automatic weather station (AWS) (Adcon Telemetry, A733 GSM/GPRS, Klosterneuburug, Austria) was installed over a well-watered grass at 1000 m to the northeast of the apple orchard to measure the solar radiation (R s ) (W•m −2 ), wind speed (u 2 ) (m•s −1 ), relative humidity (RH) (%), and air temperature (T a ) ( • C) at 15-min time intervals (Figure 1b).These variables were used as inputs to compute instantaneous (ETr i ) and daily reference evapotranspiration (ETr) using FAO Penman-Monteith equation [2].

Energy Balance Component Measurements
A tower of 5.5 m height was installed to measure the energy balance components over the dripirrigated apple orchard (Figure 1c).Latent (LEE) and sensible (HE) heat fluxes were measured with an eddy covariance (EC) system, which is composed of a three-dimensional sonic anemometer (CSAT3, Campbell Scientific Inc., Logan, UT, USA) and an infrared gas analyzer (LI-7500, LI-COR Inc., Lincoln, NE, USA).Instantaneous measurements of temperature and humidity were registered at a

Energy Balance Component Measurements
A tower of 5.5 m height was installed to measure the energy balance components over the drip-irrigated apple orchard (Figure 1c).Latent (LE E ) and sensible (H E ) heat fluxes were measured with an eddy covariance (EC) system, which is composed of a three-dimensional sonic anemometer (CSAT3, Campbell Scientific Inc., Logan, UT, USA) and an infrared gas analyzer (LI-7500, LI-COR Inc., Lincoln, NE, USA).Instantaneous measurements of temperature and humidity were registered at a frequency of 10 Hz and stored on a data logger (CR5000, Campbell Scientific Inc., Logan, UT, USA).The post processing of raw EC data included density correction [32], the sonic temperature correction [33], and coordinate rotation [34].
The net radiation (Rn E ) was measured using a four-way net radiometer (CNR1, Kipp and Zonen Inc., Delft, The Netherlands) and soil heat flux (G E ) was obtained using soil heat flux plates buried at a depth of 0.08 m (HFT3, Campbell Scientific Inc., Logan, UT, USA).In addition, soil temperature was measured by four soil thermocouples (TCAV, Campbell Scientific Inc., Logan, UT, USA) positioned at 0.02 and 0.06 m depth above each heat flux plate [35].Four plates were installed under tree rows and four other plates were located between rows.
The method proposed by [36] (R92 model) was used to estimate Z om .In this case, the wind velocities from the CSAT3 and canopy height were used to compute Z om values over the EC footprint area.Data processing, equations and theoretical basis used in this study are described by [37].

Data Quality Control and Post-Processing
To reduce the uncertainty associated with error in the LE E and H E measurements, entire days were excluded from the study when the ratios of (H E + LE E ) to (Rn E − G E ) were outside the range of 0.8 and 1.2 [7].In addition, H E and LE E values were adjusted using the Bowen ratio (β = H E /LE E ) to close the energy balance for measured fluxes as follows [38][39][40]:

Satellite Images and Processing Procedure
Twenty Landsat 7 ETM+ satellite images (path 233/row85) were downloaded using the USGS Global Visualization Viewer (GLOVIS, http://glovis.usgs.gov).All scenes used in this research (Table 1) considered only days having less than 30% cloud cover at the field scale, acquired between November and March.The images had a standard level "Level 1T" correction performed by the USGS, which provides geometric and radiometric correction using ground control points and digital elevation model.

Calibration of LAI, Z om and G i Functions
In general, Z om , LAI and G i functions in METRIC were developed for annual crops and may not apply well to heterogeneous or complex canopies such as orchards.In this experiment, the Weibull function was calibrated to estimate LAI as follows (determination coefficient (R 2 ) = 0.87): LAI = 2.42 − 1.04e −502.1•NDVI 9.32  (11) For heterogeneous canopies, the Perrier function [41] has been suggested to estimate Z om in olive orchards, vineyards and forestry plantations [6,18,23].The Perrier function is expressed as (R 2 = 0.82): where h (m) is the crop height; a is an adjustment factor describing the LAI distribution within the canopy.In this study, a was 0.06, rather than the standard values of 0.4 to 0.8 [6] due to the impacts of training system and canopy architecture of the apple orchard.Finally, G i was computed using the following recalibrated equation:

Model Validation
The validation of sub-models that compute the instantaneous energy balance components (Rn i , G i , LE i , and H i ), LAI and Z om is described using the ratio of estimated to observed values (b), mean absolute error (MAE), root-mean-square error (RMSE) and index of agreement (d).The student's t-test analysis was applied to check whether b was significantly different from unity at the 95% confidence level.Finally, hourly measurements of Rn i , G i , H β and LE β between 11:00 and 12:00 h (local time) were linearly interpolated to obtain the instantaneous values at the time of the satellite overpass (about 11:30 h).Values of RMSE, MAE and d were respectively computed as [42,43]: where n is the number of observations, y i is the observed value, and ŷi is the modeled value.Finally, the model validation was carried out using field data collected during in the 2012/2013 season, while model calibration was carried out using ground-truth data collected during 2013/2014 and 2014/2015 seasons (see Table 1).

Climatic and Crop Conditions
Predominantly clear sky days were observed during the three study periods (from November to April) with dry and hot atmospheric conditions.Maximum values for T a , vapor pressure deficit (VPD) and ETr were observed at satellite overpass time from mid-December to early January (local summer), with means of 28.7 • C (±3.9 • C), 9.1 kPa (±3.8 kPa), and 6.1 mm•day −1 (±0.9 mm•day −1 ), respectively.Cumulative ETr from September to April ranged between 829 and 904 mm with a maximum value observed during 2013/2014 season.Annual rainfall was mainly concentrated from May to September and totaled 585, 601, and 707 mm for the 2012/2013, 2013/2014 and 2014/2015 study periods, respectively, with no rainfall events near the time of satellite overpasses (Figure 2).The volumetric soil water content at rooting depth and Ψ md ranged between 0.24 and 0.33 m 3 •m −3 , and −1.01 and −1.37 MPa, respectively, indicating that the high-density apple orchard was not under water stress conditions during the three study periods [44].
Remote Sens. 2017, 9, 670 8 of 18 volumetric soil water content at rooting depth and Ψmd ranged between 0.24 and 0.33 m 3 •m −3 , and −1.01 and −1.37 MPa, respectively, indicating that the high-density apple orchard was not under water stress conditions during the three study periods [44].
The energy balance closure for hourly flux measurements is indicated in Figure 3 which shows that many observations were close to the 1:1 line.Statistical analyses indicated that the slope of the regression line through the origin was statistically different from unity at the 95% confidence level, and suggested that turbulent fluxes (HE + LEE) were less than the available energy (RnE − GE) by about 14%.According to several researchers, this result is considered appropriate to provide accurate estimates of HE + LEE, especially following adjustment using Equations ( 9) and ( 10) [45][46][47][48][49].The energy balance closure for hourly flux measurements is indicated in Figure 3 which shows that many observations were close to the 1:1 line.Statistical analyses indicated that the slope of the regression line through the origin was statistically different from unity at the 95% confidence level, and suggested that turbulent fluxes (H E + LE E ) were less than the available energy (Rn E − G E ) by about 14%.According to several researchers, this result is considered appropriate to provide accurate estimates of H E + LE E , especially following adjustment using Equations ( 9) and ( 10) [45][46][47][48][49].

Validation of the Original and Calibrated Sub-Models of LAI and Zom
In this study, LAI measured in the apple orchard ranged between 1.65 and 2.41 m 2 •m −2 while LAI estimated from original standard equations in METRIC (Equations ( 4) and ( 5 2).In addition, the t-test indicated that b was significantly equal to unity indicating that estimated and observed values of LAI were similar.In contrast, the original approaches (Equations ( 4) and ( 5)) underestimated LAI by about 59% with RMSE and MAE equal to 1.24 and 1.16 m 2 •m −2 , respectively (Table 3).In this regard, Figure 4 indicates that all points from the original standard algorithms were below the line 1:1 while those from the calibrated Weibull model were close to the 1:1 line.Some of the underestimation by the standard functions was due to the tall, narrow nature of the pruned apple trees, where NDVI and SAVI from the nadir-viewing Landsat satellite underestimated actual LAI.Carrasco-Benavides et al. (2014) [18] showed that the performance of the Weibull model depends on the shape of the canopy, fractional cover, and cover crop maintained between rows.In this study, a constant value of fraction cover (0.35 ± 0.03) was observed and the soil surface between rows was kept without weeds, especially from December to March.Major disagreements between observed and estimated values of LAI were found during mid spring (November) when the soil surface between rows was partially covered by weeds.Under this condition, values for NDVI and SAVI were overestimated because satellite images (30 m resolution) failed to represent only orchard cover.Hall et al. (2008) [50] described how the moderate spatial resolution from satellites (such as Landsat) produce mixed-pixel information for vineyards that includes both row and inter-row space information.Agreement between observed and estimated LAI improved after early December when weeds were cut and the soil surface between rows became dry.Finally, the wetted area of the drip-irrigated orchard was only about 5.0% of the total area and was located below the tree canopy.

Validation of the Original and Calibrated Sub-Models of LAI and Z om
In this study, LAI measured in the apple orchard ranged between 1.65 and 2.41 m 2 •m −2 while LAI estimated from original standard equations in METRIC (Equations ( 4) and ( 5)) ranged between 0.  2).In addition, the t-test indicated that b was significantly equal to unity indicating that estimated and observed values of LAI were similar.In contrast, the original approaches (Equations ( 4) and ( 5)) underestimated LAI by about 59% with RMSE and MAE equal to 1.24 and 1.16 m 2 •m −2 , respectively (Table 3).In this regard, Figure 4 indicates that all points from the original standard algorithms were below the line 1:1 while those from the calibrated Weibull model were close to the 1:1 line.Some of the underestimation by the standard functions was due to the tall, narrow nature of the pruned apple trees, where NDVI and SAVI from the nadir-viewing Landsat satellite underestimated actual LAI.Carrasco-Benavides et al. (2014) [18] showed that the performance of the Weibull model depends on the shape of the canopy, fractional cover, and cover crop maintained between rows.In this study, a constant value of fraction cover (0.35 ± 0.03) was observed and the soil surface between rows was kept without weeds, especially from December to March.Major disagreements between observed and estimated values of LAI were found during mid spring (November) when the soil surface between rows was partially covered by weeds.Under this condition, values for NDVI and SAVI were overestimated because satellite images (30 m resolution) failed to represent only orchard cover.Hall et al. (2008) [50] described how the moderate spatial resolution from satellites (such as Landsat) produce mixed-pixel information for vineyards that includes both row and inter-row space information.Agreement between observed and estimated LAI improved after early December when weeds were cut and the soil surface between rows became dry.Finally, the wetted area of the drip-irrigated orchard was only about 5.0% of the total area and was located below the tree canopy.This study indicates that the calibrated Perrier equation (Equation ( 12)) using a = 0.06 was able to predict Z om with RMSE and MAE both equal to 0.03 m (Table 2).In addition, b was significantly equal to unity, indicating that estimated and observed values for Z om were similar.Furthermore, the original standardized equation (Equation ( 8)) underestimated Z om by about 83% with RMSE = 0.22 m and MAE = 0.21 m (Table 3).Figure 5 shows that the estimates by the Perrier equation are close to the 1:1 line when a = 0.06, but those for the standard model setting (a = 0.6) were below the straight line.These results are similar to those found by [18] who indicated that Z om was underestimated by about 81% when the original function setting was applied to a drip-irrigated vineyard.For ETa estimation from a non-irrigated olive orchard, Santos et al. (2012) [23] found that the implementation of the Perrier function reduced the values of RMSE from 1.12 mm•day −1 to 0.25 mm•day −1 .In addition, Pôças, et al. (2014) [22] observed that the Perrier equation applied on a dense hedgerow olive orchard produced similar results to those found in this study where Z om values varied with LAI and canopy height.These authors also indicated that original equation settings utilized in METRIC for computing Z om (Equation ( 8)) tend to produce an underestimation on H i in the orchard.
This study indicates that the calibrated Perrier equation (Equation ( 12)) using a = 0.06 was able to predict Zom with RMSE and MAE both equal to 0.03 m (Table 2).In addition, b was significantly equal to unity, indicating that estimated and observed values for Zom were similar.Furthermore, the original standardized equation (Equation ( 8)) underestimated Zom by about 83% with RMSE = 0.22 m and MAE = 0.21 m (Table 3).Figure 5 shows that the estimates by the Perrier equation are close to the 1:1 line when a = 0.06, but those for the standard model setting (a = 0.6) were below the straight line.These results are similar to those found by [18] who indicated that Zom was underestimated by about 81% when the original function setting was applied to a drip-irrigated vineyard.For ETa estimation from a non-irrigated olive orchard, Santos et al. (2012) [23] found that the implementation of the Perrier function reduced the values of RMSE from 1.12 mm•day −1 to 0.25 mm•day −1 .In addition, Pôças, et al. (2014) [22] observed that the Perrier equation applied on a dense hedgerow olive orchard produced similar results to those found in this study where Zom values varied with LAI and canopy height.These authors also indicated that original equation settings utilized in METRIC for computing Zom (Equation ( 8)) tend to produce an underestimation on Hi in the orchard.

Validation of Orchard Energy Balance Fluxes
The comparisons between observed and estimated values for LEi, Rni, Hi and Gi at the time of satellite overpasses (11:30-11:40 h local time) are depicted in Figure 6a,b when using the original and calibrated sub-models to compute Zom, LAI and Gi, respectively.For net radiation (triangles), these figures show that all points were close to the 1:1 line, with RMSE and MAE ranging 15-20 and 15-18 W•m −2 , respectively (Tables 2 and 3).For both original and calibrated sub-models, the t-test indicated that b was not different from unity (one) indicating that observed and estimated values of Rni were similar.In this study, a variation of ±50% in the estimation of LAI produced differences of ±3% in the calculation of broad-band surface emissivity (εo = 0.95 + 0.01 LAI) which is used to estimate the outgoing long-wave radiation (Equation ( 2)).This sensitivity analysis explains why the calibrated version of Weibull approach did not improve the calculation of Rni.
The comparison between observed and estimated Gi indicates that the calibration of c1 and c2 coefficients reduced the error from 12% to 3%.Values for RMSE and MAE were 33 and 23 W•m −2 for the original model while those were 16 and 14 W•m −2 for the calibrated function, respectively.In addition, Figure 6a,b shows that the points were grouped around the 1:1 line.Carrasco-Benavides et al. (2014) [18] and Samani and Bawazir (2015) [51] suggested that error in the estimation of Gi can be associated with the different measurement scales between satellite data and ground-based observations and sun angle.

Validation of Orchard Energy Balance Fluxes
The comparisons between observed and estimated values for LE i , Rn i , H i and G i at the time of satellite overpasses (11:30-11:40 h local time) are depicted in Figure 6a,b when using the original and calibrated sub-models to compute Z om , LAI and G i , respectively.For net radiation (triangles), these figures show that all points were close to the 1:1 line, with RMSE and MAE ranging 15-20 and 15-18 W•m −2 , respectively (Tables 2 and 3).For both original and calibrated sub-models, the t-test indicated that b was not different from unity (one) indicating that observed and estimated values of Rn i were similar.In this study, a variation of ±50% in the estimation of LAI produced differences of ±3% in the calculation of broad-band surface emissivity (ε o = 0.95 + 0.01 LAI) which is used to estimate the outgoing long-wave radiation (Equation ( 2)).This sensitivity analysis explains why the calibrated version of Weibull approach did not improve the calculation of Rn i .
The comparison between observed and estimated G i indicates that the calibration of c 1 and c 2 coefficients reduced the error from 12% to 3%.Values for RMSE and MAE were 33 and 23 W•m −2 for the original model while those were 16 and 14 W•m −2 for the calibrated function, respectively.In addition, Figure 6a,b shows that the points were grouped around the 1:1 line.Carrasco-Benavides et al. (2014) [18] and Samani and Bawazir (2015) [51] suggested that error in the estimation of G i can be associated with the different measurement scales between satellite data and ground-based observations and sun angle.Results indicated that METRIC using the original, standard functions underestimated Hi by about 29% with RMSE and MAE ranging 80-95 W•m −2 (Table 3).In this case, Figure 6a shows that points were below the straight line with an important dispersion.When using the calibrated functions for Zom, METRIC overestimated Hi by about 5% with RMSE = 33 W•m −2 and MAE = 26 W•m −2 (Table 2).For these results, Figure 6b shows that points were close to the 1:1 line.In addition, underlying error found in this study may be associated with the selection of the hot pixels which were difficult to define accurately for some satellite scenes.Thus, the internal self-calibration of METRIC (called CIMEC) plays a critical role in achieving accurate values for LEi because it absorbs many errors and biases generated during the computation of Rni and Gi when Hi is estimated [16].
Finally, the statistical analyses indicate that METRIC using the original standard functions was able to simulate LEi with RMSE and MAE of 95 and 81 W•m −2 , respectively (Table 3).The slope b was statistically greater than one indicating that LEi was overestimated by about 26%.When using the calibrated functions, METRIC overestimated LEi over the high-density apple orchard by about 16% with a RMSE = 30 W•m −2 , MAE = 20 W•m −2 and d = and 0.87 (Table 3).Figure 6 (rhombs) shows that points were closer to the 1:1 when using METRIC with the recalibrated functions for Zom, LAI and Gi.These results suggest that the calibration of LAI, Zom, and Gi improved the performance of METRIC to estimate LEi over the drip-irrigated apple orchard.Allen et al. (2011) [3] indicated that one of the main disadvantages of METRIC is that the accuracy in estimating LEi for individual pixels depends on errors in the calculations for Rni, Hi, and Gi.In general, the smaller the Zom computed by the original sub-model of METRIC, the lower the values for Hi, and hence estimations of LEi are typically overestimated.This was the case for the apple orchards.This behavior is particularly important considering that Zom values are typically underestimated in fruit orchards and vineyards using the standard sub-model parameters (Equation ( 8)) and consequently the daily actual evapotranspiration is overestimated [23].
Results observed in this study are similar to those found in the literature [18] that indicated that METRIC estimated LEi with errors of 15% in a drip-irrigated vineyard.González-Dugo et al. (2012) [52] reported an error of 18% in the estimation of LEi over vine grapes.In an orange orchard, LEi fluxes computed by a one-source remote sensing model were similar to those measured by an Eddy Covariance system at the field scale, with RMSE equal to 25 W•m −2 [15].Usually, the overall mean bias in estimated LEi and ETa varies between 2% and 25% using one-source models for different fruit crops [53].These results are considered acceptable by many researchers indicating potential to evaluate different irrigation management with remote sensing techniques and to improve irrigation efficiency at different field scales [6,19,20,54].Results indicated that METRIC using the original, standard functions underestimated H i by about 29% with RMSE and MAE ranging 80-95 W•m −2 (Table 3).In this case, Figure 6a shows that points were below the straight line with an important dispersion.When using the calibrated functions for Z om , METRIC overestimated H i by about 5% with RMSE = 33 W•m −2 and MAE = 26 W•m −2 (Table 2).For these results, Figure 6b shows that points were close to the 1:1 line.In addition, underlying error found in this study may be associated with the selection of the hot pixels which were difficult to define accurately for some satellite scenes.Thus, the internal self-calibration of METRIC (called CIMEC) plays a critical role in achieving accurate values for LE i because it absorbs many errors and biases generated during the computation of Rn i and G i when H i is estimated [16].
Finally, the statistical analyses indicate that METRIC using the original standard functions was able to simulate LE i with RMSE and MAE of 95 and 81 W•m −2 , respectively (Table 3).The slope b was statistically greater than one indicating that LE i was overestimated by about 26%.When using the calibrated functions, METRIC overestimated LE i over the high-density apple orchard by about 16% with a RMSE = 30 W•m −2 , MAE = 20 W•m −2 and d = and 0.87 (Table 3).Figure 6 (rhombs) shows that points were closer to the 1:1 when using METRIC with the recalibrated functions for Z om , LAI and G i .These results suggest that the calibration of LAI, Z om , and G i improved the performance of METRIC to estimate LE i over the drip-irrigated apple orchard.Allen et al. (2011) [3] indicated that one of the main disadvantages of METRIC is that the accuracy in estimating LE i for individual pixels depends on errors in the calculations for Rn i , H i , and G i .In general, the smaller the Z om computed by the original sub-model of METRIC, the lower the values for H i , and hence estimations of LE i are typically overestimated.This was the case for the apple orchards.This behavior is particularly important considering that Z om values are typically underestimated in fruit orchards and vineyards using the standard sub-model parameters (Equation ( 8)) and consequently the daily actual evapotranspiration is overestimated [23].
Results observed in this study are similar to those found in the literature [18] that indicated that METRIC estimated LE i with errors of 15% in a drip-irrigated vineyard.González-Dugo et al. (2012) [52] reported an error of 18% in the estimation of LE i over vine grapes.In an orange orchard, LE i fluxes computed by a one-source remote sensing model were similar to those measured by an Eddy Covariance system at the field scale, with RMSE equal to 25 W•m −2 [15].Usually, the overall mean bias in estimated LE i and ETa varies between 2% and 25% using one-source models for different fruit crops [53].These results are considered acceptable by many researchers indicating potential to evaluate different irrigation management with remote sensing techniques and to improve irrigation efficiency at different field scales [6,19,20,54].

Distribution of Energy Balance into the Apple Orchard
Maps showing the spatial and temporal variability of the energy balance components over a drip-irrigated apple orchard are illustrated in Figure 7.For the 2012/2013 growing season, the evolution of instantaneous net radiation (Rn i ) (Figure 7a-g), soil heat flux (G i ) (Figure 7h-o), sensible heat flux (H i ) (Figure 7p-v), and latent heat flux (LE i ) (Figure 7w-cc) were computed using METRIC with the calibrated sub-models of Z om , LAI and G i .In the study area, mean values for Rn i ranged between 494 and 636 W•m −2 with minimum and maximum values observed on DOY 78 and 316, respectively.The spatial variability in Rn i was quite homogeneous over the apple orchard with a coefficient of variation less than 2% during the study period.Maximum and minimum values for G i were found on DOY 364 (dry and hot atmospheric conditions) and DOY 78 (beginning of autumn in the southern hemisphere) with average values of 110 (±6) and 82 (±10) W•m −2 , respectively.For G i , the spatial distribution exhibited a coefficient of variation less than 11% which agreed with changes observed in the scattering during the season.When the fruit reached half of their size (DOY 364), H i exhibited the highest spatial variability observed over the orchard with mean values of 273 (±44) W•m −2 .In contrast, at the beginning of ripening (DOY46), the lowest spatial variability was detected with mean values of 259 (±3) W•m −2 and with a coefficient of variation of 3%.Because of the energy balance budget, the spatial variability of LE i exhibited a coefficient of variation of 4% and 26% for DOY 364 and DOY 46, respectively with mean values of 200 (±52) and 287 (±10) W•m −2 in each case.At the beginning of the season (DOY 316), the highest values for LE i were observed in the south side of the plot where a grass cover between rows was growing due to winter rainfall and rising temperatures (mid spring).This condition coincided with maximum differences observed between measured and modeled energy balance fluxes for this day, which led into an overestimation in LE i by the model.During the pre-harvest period (DOY 78), irrigation water supply was partially suspended in the orchard to predispose fruits and trees for harvest a month later.This issue can be observed as a relatively homogeneous spatial distribution of energy balance components through the orchard, according to the remote sensing analysis.
Application of RSEB models like METRIC to compute orchard water requirements has recently been applied with promising results [22,23].However, the main limitation of METRIC is that the estimation of LE i for each pixel is only as accurate as the estimates for Rn i , G i and H i [6].In this study, the main errors were associated with the selection of hot pixels within the scene which were very difficult to obtain for some satellite images.This limitation may cause significant errors in the computation of LE i by the model, especially in heterogeneous crops, like grapes and apples, where H could generate a big impact in the orchard energy balance (Equation ( 1)).In this regard, López-Olivarí et al. (2016) [7] indicated that the canopy architecture (tree density, canopy size, LAI, and fc) and training system can have an important effect in the partitioning of Rn and ETa inside an orchard and in the distribution of the energy balance components over the land surface.In addition, Testi et al. (2004) [48] indicated that the exchange of energy fluxes (Rn i , G i , H i and LE i ) between soil surface and tree canopy is strongly affected by the orchard architecture and its shading, and the sensible heat flux generated at the soil surface can be a major contributor to the energy balance for sparse woody canopies playing a key role in tree transpiration and stomatal closure [55].However, this issue is not possible to demonstrate accurately in METRIC due to the limitations of one-source models to compute the surface energy balance fluxes using satellite-based blended data of soil and canopy at the 30 × 30 m 2 pixel size.In addition, main limitations of current satellite platforms for practical application of site-specific irrigation management, especially in heterogeneous canopies such as orchards and vineyards, are the lack of fine spatial resolution and real-time data at the field and sub-field scales [56].

Conclusions
METRIC was originally developed to estimate the land surface energy balance and ETa over relatively homogeneous crop canopies, mostly associated with annual crops having full ground cover.However, over sparse woody canopies, standard functions in METRIC tend to produce underestimations for LAI, Z om and G i fluxes.This issue generates lower values for H i and hence high LE i production via the energy balance budget.After adopting calibrated algorithms for LAI, Z om , and G i , the intermediate components of land surface energy balance were adjusted, and error in LAI, G i and Z om values over the apple orchard was substantially reduced.
Results suggest that METRIC can be used as an effective tool to simulate the spatial and temporal variability of energy balance components with relative accuracy for a high-density apple orchard, with error within the range of biases observed by other authors who have studied tree crops.However, it is necessary to conduct more research in additional orchards and with different types of fruit trees in order to analyze the behavior of intermediate algorithms of the model, especially in counties such as Chile were a single Landsat scene (180 km 2 ) contains three different microclimatic conditions: The Andes Mountains, the central valley, and the coastal range.Future research should also consider the implementation of the two-source or three-source energy balance models to analyze the behavior of energy balance components over orchards.

Figure 1 .
Figure 1.Location of the study area in the Pelarco valley, Maule Region (Central Chile), with the identification of: (a) the super intensive apple orchard; (b) the automatic weather station; and (c) the Eddy Covariance system.

Figure 1 .
Figure 1.Location of the study area in the Pelarco valley, Maule Region (Central Chile), with the identification of: (a) the super intensive apple orchard; (b) the automatic weather station; and (c) the Eddy Covariance system.

Figure 2 .
Figure 2. Daily values of reference evapotranspiration (ETr) and precipitation events (vertical solid lines) over the high-density apple orchard during the study seasons (Pelarco valley, Maule Region, Chile).Satellite overpass dates are indicated as vertical dashed lines.

Figure 2 .
Figure 2. Daily values of reference evapotranspiration (ETr) and precipitation events (vertical solid lines) over the high-density apple orchard during the study seasons (Pelarco valley, Maule Region, Chile).Satellite overpass dates are indicated as vertical dashed lines.

Figure 4 .
Figure 4. Comparison between observed and estimated values of Leaf Area Index (LAI) of dripirrigated apple orchard (Pelarco Valley, Maule Region, Chile): (a) original short crop function; and (b) calibrated by Weibull function for summer periods.
mean square error; MAE = mean absolute error; d = index of agreement; b = slope between estimated/observed T = true hypothesis (b = 1); F = false hypothesis (b ≠ 1).Subscripts: "i" indicated instantaneous values computed by METRIC model, "E" is measured by the Eddy covariance system, "β" denotes fluxes corrected by Bowen ratio, "w" is computed by Weibull calibrated function and "C" denotes calibrated coefficients, respectively.
mean square error; MAE = mean absolute error; d = index of agreement; b = slope between estimated/observed T = true hypothesis (b = 1); F = false hypothesis (b ≠ 1).Subscripts: "i" indicated instantaneous values computed by METRIC model "M" is computed by original METRIC model, "E" is measured by the Eddy covariance system and "β" denotes fluxes corrected by Bowen ratio, respectively.

Figure 4 .
Figure 4. Comparison between observed and estimated values of Leaf Area Index (LAI) of drip-irrigated apple orchard (Pelarco Valley, Maule Region, Chile): (a) original short crop function; and (b) calibrated by Weibull function for summer periods.

Figure 5 .
Figure 5.Comparison between observed and estimated values for Aerodynamic Roughness length (Zom) for a drip-irrigated apple orchard (Pelarco Valley, Maule Region, Chile): (a) original short crop function; and (b) adapted Perrier function.

Figure 5 .
Figure 5.Comparison between observed and estimated values for Aerodynamic Roughness length (Z om ) for a drip-irrigated apple orchard (Pelarco Valley, Maule Region, Chile): (a) original short crop function; and (b) adapted Perrier function.

Figure 6 .
Figure 6.Comparisons between measured and estimated values for instantaneous latent heat flux (LEi), net radiation (Rni), sensible heat flux (Hi) and soil heat flux (Gi) at the time of satellite overpass (11:30-11:40 h local time) during the 2012/2013 growing season for the apple orchard (Pelarco Valley, Maule Region, Chile): (a) original sub-models; and (b) calibrated sub-models.

Figure 6 .
Figure 6.Comparisons between measured and estimated values for instantaneous latent heat flux (LE i ), net radiation (Rn i ), sensible heat flux (H i ) and soil heat flux (G i ) at the time of satellite overpass (11:30-11:40 h local time) during the 2012/2013 growing season for the apple orchard (Pelarco Valley, Maule Region, Chile): (a) original sub-models; and (b) calibrated sub-models.

Figure 7 .
Figure 7. Spatial and temporal variability of: instantaneous net radiation (Rni) (a-g) soil heat flux (Gi) (h-o); sensible heat flux (Hi) (p-v); and latent heat flux (LEi) (w-cc) computed using calibrated algorithms of METRIC model over the study site at the time of satellite flybys (11:30-11:40 h local time), during the 2012/2013 growing season.

Figure 7 .
Figure 7. Spatial and temporal variability of: instantaneous net radiation (Rn i ) (a-g) soil heat flux (G i ) (h-o); sensible heat flux (H i ) (p-v); and latent heat flux (LE i ) (w-cc) computed using calibrated algorithms of METRIC model over the study site at the time of satellite flybys (11:30-11:40 h local time), during the 2012/2013 growing season.

Table 1 .
Summary of Landsat 7 ETM+ satellite imagery used to estimate land surface energy balance fluxes over the apple orchard (path 233/row 85).

Table 2 .
Statistical analysis results from recalibrated sub-models to estimate leaf area index (LAI), momentum roughness length (Zom), net radiation (Rni), soil heat flux (Gi), sensible heat flux (Hi), and latent heat flux (LEi) by the METRIC model over the high-density apple orchard for the 2012/2013 growing season, Pelarco Valley, Chile.

Table 3 .
Statistical analysis results from the original sub-models to estimate leaf area index (LAI), momentum roughness length (Zom), net radiation (Rni), soil heat flux (Gi), sensible heat flux (Hi), and latent heat flux (LEi) by the METRIC model over the high-density apple orchard for the 2012/2013 growing season, Pelarco Valley, Chile.

Table 2 .
Statistical analysis results from recalibrated sub-models to estimate leaf area index (LAI), momentum roughness length (Z om ), net radiation (Rn i ), soil heat flux (G i ), sensible heat flux (H i ), and latent heat flux (LE i ) by the METRIC model over the high-density apple orchard for the 2012/2013 growing season, Pelarco Valley, Chile.Using Recalibrated LAI, Z

om and G i Algorithms of METRIC
RMSE = root mean square error; MAE = mean absolute error; d = index of agreement; b = slope between estimated/observed T = true hypothesis (b = 1); F = false hypothesis (b = 1).Subscripts: " i " indicated instantaneous values computed by METRIC model, " E " is measured by the Eddy covariance system, " β " denotes fluxes corrected by Bowen ratio, " w " is computed by Weibull calibrated function and " C " denotes calibrated coefficients, respectively.

Table 3 .
Statistical analysis results from the original sub-models to estimate leaf area index (LAI), momentum roughness length (Z om ), net radiation (Rn i ), soil heat flux (G i ), sensible heat flux (H i ), and latent heat flux (LE i ) by the METRIC model over the high-density apple orchard for the 2012/2013 growing season, Pelarco Valley, Chile.
RMSE = root mean square error; MAE = mean absolute error; d = index of agreement; b = slope between estimated/observed T = true hypothesis (b = 1); F = false hypothesis (b = 1).Subscripts: " i " indicated instantaneous values computed by METRIC model " M " is computed by original METRIC model, " E " is measured by the Eddy covariance system and " β " denotes fluxes corrected by Bowen ratio, respectively.