Evaluation of transpiration in different almond production systems with two-source energy balance models from UAV thermal and multispectral imagery

Abstract

A growing number of intensive irrigated production systems of the almond crop have been established in recent years. However, there is little information regarding the crop water requirements. Remote sensing-based models such as the two-source energy balance (TSEB) have proven to be reliable ways to accurately estimate actual crop evapotranspiration. However, few efforts have been made to validate the transpiration with sap flow measurements in woody row crops with different production systems and water status. In this study, the TSEB Priestley-Taylor (TSEB-PT) and contextual approach (TSEB-2T) models were assessed to estimate canopy transpiration. In addition, the effect of applying a basic clumping index for heterogeneous randomly placed clumped canopies and a rectangular hedgerow clumping index on the TSEB transpiration estimation was assessed. The TSEB inputs were obtained from high resolution multispectral and thermal imagery using an unmanned aerial vehicle. The leaf area index (LAI), stem water potential (Ψ_stem) and fractional intercepted photosynthetically active radiation (fIPAR) were also measured. Significant differences were observed in transpiration between production systems and irrigation treatments. The combined use of the TSEB-2T with the C&N-R transmittance model gave the best transpiration estimations for all production systems and irrigation treatments. The use of in situ PAR transmittance in the TSEB-2T model significantly improved the root mean squared error. Thus, the better agreement observed with the TSEB when using the C&N-R model and in situ PAR transmittance highlights the importance of improving radiative transfer models for shortwave canopy transmittance, especially in woody row crops.

Introduction

Almond production has increased worldwide in the last decade (FAOSTAT 2022). Currently, more than 2 million ha of orchards are cultivated worldwide (FAOSTAT 2022). In Spain, almond has traditionally been planted in rainfed marginal areas with poor and shallow soils. However, the establishment of modern orchards in new irrigated areas has led to significant progress in the production techniques of almond tree cultivation. Among them, intensive plantations, characterized by smaller spacing distances and more planar canopies, are increasingly common. The area of irrigated almond cultivation in Spain has increased from 4.7% in 2005 to 22.1% in 2021. Moreover, new intensive orchards amounted to over 140,000 ha in 2020 (MAPA 2021). This trend in orchard intensification from 3D canopy architectures to modern high density, simple/planar designs is happening in other Prunus species as well as almond (Iglesias and Echeverria 2022).

The shift to a more intensified agriculture coincides in time with a context of water scarcity in which water resources are already limited in, for example, the Mediterranean countries (Tramblay et al. 2020; Moldero et al. 2021, 2022; Soares and Lima 2022). Therefore, proper water management of irrigated crops will be essential to ensure successful long-term agricultural activity (Garcia-Tejero et al. 2014). In this context, the adoption of regulated deficit irrigation (RDI) strategies plays an important role in contributing to reducing water consumption without significantly impacting yield (Girona et al. 2005; Egea et al. 2010). Numerous studies have quantified crop water requirements and analyzed the impact of RDI strategies in almond (Girona et al. 2005; Fereres and Soriano 2007; Espadafor et al. 2017; López-López et al. 2018a, b; Moldero et al. 2021, 2022). However, these studies have mostly considered open-vase training systems with large spacing distances. Thus, little is known about actual crop water usage of the new and more intensive almond production systems with planar canopies or their response to the adoption of RDI strategies. Transpiration or crop water use is mostly driven by the atmospheric saturation deficit, the amount of solar radiation intercepted by the canopy, and regulated by stomatal and aerodynamic conductance (Ayars et al. 2003; Girona et al. 2011; Espadafor et al. 2015). The balance between the last two components depends on the degree of coupling of leaves to atmosphere (Jarvis 1985). The implantation of novel training systems at/or different planting densities also has a significant impact on light interception and therefore on transpiration. For instance, Casanova-Gascón et al. (2019) reported that almonds in an open-center system resulted in higher light intercepted than a superintensive system. Iglesias and Echeverria (2022) also showed that although the intercepted light of planar canopies was compensated for by a higher tree density, less light was intercepted than with other 3D canopy systems like the Y-trellis, open vase, transversal ypsilon or double Y. Despite the differences encountered in light interception, to the best of our knowledge no evapotranspiration comparisons have been made between different trellis systems.

A number of methodologies have been developed to measure or estimate crop evapotranspiration or its components (transpiration and soil evaporation) across a spectrum of spatial scales and crops. Examples of measurement techniques include weighing lysimeters (Girona et al. 2011; López-Urrea et al. 2012), eddy covariance (Bellvert et al. 2018; Drechsler et al. 2022; Knipper et al. 2023), soil water balance (López-López et al. 2018a; Moldero et al. 2021) and sap flows (Espadafor et al. 2015; Mancha et al. 2021). Sap flow sensors offer significant advantages, measuring the transpiration of each plant in a simple, continuous, automated way, and with high temporal resolution (Smith and Allen 1996; Forster 2017; Fernandez et al. 2001). Among the sap flow measuring methods available, the Compensation Heat-Pulse (CHF) has been suggested as a tool for detecting water stress and for irrigation scheduling purposes (Fernandez et al. 2001; Alarcón et al. 2005). Nonetheless, it is advisable to calibrate measurements to accurately assess transpiration for the entire plant, given the azimuthal variability in sap velocity within the trunk (López-Bernal et al. 2010; Forster 2017; Noun et al. 2022).

In recent decades, remote sensing surface energy balance models have also been used to estimate actual crop evapotranspiration (ETa) in a wide variety of environments and ecosystems (Shuttleworth and Wallace 1985; Bastiaanssen et al. 1998; Drexler et al. 2004; Overgaard et al. 2006; Allen et al. 2007; Timmermans et al. 2007; Kalma et al. 2008; Kustas and Anderson 2009; Gómez-Candón et al. 2021). Several energy balance model schemes have been developed; these can be classified as either one-source or two-source schemes (Bastiaanssen et al. 1998; Allen et al. 2007; Kalma et al. 2008; Kustas and Anderson 2009). In the one-source approach, the land surface is treated as if it were one large leaf with a single uniform layer, with no distinction being made between sinks related to vegetation and soil. The Land Surface Temperature (LST) in the scene is used to define the upper (LE ~ Rn -G) and lower (LE ~ 0) limits of evapotranspiration (Bastiaanssen et al. 1998; Allen et al. 2007). In contrast, the two-source approach parameterizes the biophysical processes to separately estimate plant transpiration and evaporation from non-vegetated surfaces (Norman et al. 1995; Colaizzi et al. 2014; Nieto et al. 2019). The two-source scheme emerged as a strategy to address crucial factors affecting the relationship between aerodynamic and radiometric temperature, which had yielded unsatisfactory outcomes in cases of partial canopy cover and for heterogeneous landscapes when one-source modelling was used (Norman et al. 1995). Although the new one-source formulations redefine the radiometric-aerodynamic relationship, they still need similar inputs to the two-source scheme despite providing less output information, although comparable performance (Kustas and Anderson 2009; Peddinti and Kisekka 2022). Additionally, the transpiration output from the two-source scheme offers the advantage of permitting the direct determination of plant water use and canopy stress. The two-source energy balance (TSEB) has demonstrated robustness in estimating plant transpiration and evaporation fluxes for various surface conditions and many different landscapes (Kustas and Anderson 2009; Kustas et al. 2019; Gómez-Candón et al. 2021; Gao et al. 2023; Knipper et al. 2023).

The TSEB partitions the surface energy fluxes between nominal soil and canopy sources using estimates of soil (T_s) and canopy temperature (T_c). Since direct measurements of T_c and T_s cannot be directly retrieved from coarse resolution satellite images, the Priestley and Taylor (1972) formulation has been proposed to derive T_c and T_s separately from the radiometric temperature (T_rad). T_c and T_s can also be obtained separately using unmanned aerial vehicles (UAVs) and very high resolution images (Nieto et al. 2019). Interest in the use of UAVs for monitoring crop evapotranspiration or water status has increased in recent years (López-Olivari et al. 2016; Nieto et al. 2019; Niu et al. 2020; Peddinti and Kisekka 2022; Gao et al. 2023; Ramírez-Cuesta et al. 2023). Nieto et al. (2019) showed in grapevine that using the TSEB with a simple contextual algorithm to derive soil and canopy temperatures separately (TSEB-2T) yielded the closest agreement with flux tower measurements. Bellvert et al. (2018) also reported a better performance of the TSEB-2T approach in estimating vine transpiration in comparison to the Priestley-Taylor approach. More recently, Gao et al. (2023) proposed a new method for temperature partitioning based on a quantile technique separation (QTS) and high-resolution information which, coupled with the TSEB model, improved the sensible heat flux (H) estimation by 61% in comparison to the Priestley-Taylor approach.

Besides T_c and T_s, the fractional cover (f_c), canopy height (h_c), canopy width (w_c) and leaf area index (LAI) are inputs that are also required with the TSEB modelling scheme. These parameters are used to estimate the shortwave canopy transmittance and reflectance of vegetated surfaces (Campbell and Norman 1998). However, there is some uncertainty in the estimation of canopy transmittance, especially in hedge row crops with partial canopy cover. In those cases, radiative transfer models need to account for the amount of radiation that will be directly transmitted through the inter-row space as well as the radiation transmitted through canopy gaps and through the canopy leaves (Parry et al. 2019). In this regard, there have been refinements suggested to algorithms of TSEB for row crops related to radiation partitioning (Colaizzi et al. 2012; Parry et al. 2019; Nieto et al. 2019). One such suggestion is based on using the Campbell and Norman (C&N) radiation transfer model and in which the clumping index is derived considering a geometrical model with a rectangular canopy shape (C&N-R) (Parry et al. 2019; Colaizzi et al. 2012).

Determining the optimal TSEB model framework for ET partitioning in row crops remains a challenge, as is also verification of TSEB-estimated soil evaporation and canopy transpiration in trees with different canopy structures and water status. Current methodologies used to validate total ET are mostly based on eddy covariance flux towers, but evaporation and transpiration cannot be directly measured (Gao et al. 2023). On the other hand, sap flows, once calibrated, can be used as a sensor to validate actual estimates of transpiration. Therefore, the main objectives of this research are: (1) to validate and compare the estimates of almond crop transpiration made for three different production systems and irrigation treatments using sap flow measurements; and (2) to evaluate the effect of applying the shortwave transmittance model C&N-R in the TSEB using the Priestley-Taylor (PT) and contextual (2T) schemes.

Materials and methods

Trial location and design

This study was carried out during 2021 in an almond orchard located at the experimental station of the IRTA (Institute of Agrifood Research and Technology) in Les Borges Blanques, Spain (41°30′31.89’’N; 0°51′10.70’’E, 323 m elevation, Fig. 1). The study site has a Mediterranean climate and the annual precipitation and evapotranspiration recorded in 2021 were 288 and 1054 mm, respectively. During the flight campaign conducted on March 24, May 19, June 3 and 29, and July 29, 2021, the hourly air temperatures recorded during flight time were (T_a) 19.05 °C, 20.54, 28.9, 26.9 and 31.34 °C, respectively, and the hourly vapor pressure deficits (VPD) were 1.33 kPa, 1.45 kPa, 2.24 kPa, 2.33 kPa, and 2.44 kPa. The almond orchard was planted in June of 2009 under three different production systems: open vase with a spacing distance of 5.5 × 3.5 m, central axis with 5 × 3 m, and hedgerow with 4.5 × 3 m (Fig. 1). The almond scion was “Marinada” on an INRA GF-677 rootstock. Trees were managed according to commercial practices in the region. The soil texture was a clay loam and its depth ranged from 1.6 to 2 m, with field capacity and wilting point values of 27.16% and 14.32%, respectively.

Fig. 1

Location of the almond field in Les Borges Blanques (Lleida, Spain) and experimental design of the field, showing in different colours the three production systems and the three irrigation treatments (color figure online)

Trees were irrigated daily using a drip irrigation system. The open vase system was irrigated with two laterals placed on each side of the tree at 40 cm. The central axis and hedgerow systems had a lateral pipe along the row line. The open vase trees had a dripper each 70 cm with a water discharge of 2.2 L h⁻¹, while in the central axis and hedgerow systems drippers were located each 60 cm with a water discharge of 3.8 L h⁻¹ per dripper. Each production system was subjected to three different irrigation treatments: (i) Full-irrigation control, where irrigation aimed to meet ET requirements (100% ET_C) throughout the growing season, (ii) mild stress, irrigated at 50% ET_C throughout the growing season; and (iii) severe stress, irrigated at 20% ET_C throughout the growing season. Weekly irrigation was scheduled following the method described by Allen et al. (1998) which seeks to replace crop evapotranspiration (ET_C) as follows: ET_C = (ETo x Kc)—effective rainfall. ETo and Kc represent the reference evapotranspiration and crop coefficient, respectively. Effective rainfall was estimated as half of the rainfall for a single event-day with more than 10 mm of precipitation, and otherwise was considered to be zero (Olivo et al. 2009). The Penman–Monteith method was used to determine ETo (Allen et al. 1998) and the Kc values were derived from Goldhamer and Girona (2012) at different phenological stages: Kc₁ = 0.70 (April), Kc₂ = 0.95 (May), Kc₃ = 1.09 (June), Kc₄ = 1.15 (July), Kc₅ = 1.17 (August), and Kc₆ = 1.12 (September). ETo was collected from a weather station belonging to Catalonia’s official network of meteorological stations (SMC, www.ruralcat.net/web/guest/agrometeo), which was located 500 m from the study site. The amount of water applied to each treatment was measured with digital water meters (CZ2000-3M, Contazara, Zaragoza, Spain).

Field measurements

Sap flow measurement

Transpiration was estimated in situ using a sap flow measurement system based on the compensation heat pulse (CHP) method combined with the calibrated average gradient technique. The system was developed by the IAS-CSIC laboratory and corresponds to a 4.8-W stainless steel heater of 2 mm diameter and two temperature sensors located 10 and 5 mm downstream and upstream of the heater. The sap flux densities across the trunk radius are calculated on the basis of the heat pulse velocities at 5 and 15 mm below the cambium. Each temperature sensor has two embedded type E (chromel constantan wire) thermocouple junctions spaced 10 mm along the needle. For further specifications see Villalobos et al. (2009). Sap flow data were recorded every 15 min and registered in a CR1000 datalogger (Campbell Scientifc Inc., Logan, UT, USA).

In each production system, two trees of the fully irrigated and severe stress treatments and one tree of the mild stress treatment were monitored with sap flow sensors installed 0.5 m above the ground. Measurements of sap flow transpiration (T_sf) were corrected for wounding and azimuthal effects (López-Bernal et al. 2010). For each tree, sensor transpiration was corrected for actual transpiration at the beginning of the growing season using a correction coefficient. The correction coefficients were derived from in situ measurements of transpiration obtained via a water balance method using Eq. 1 and based on the strong correlation between transpiration estimated via water balance methods and all-season sap flow measurements, as described by López-López et al. (2018a).

$$T_{wb} = P + I_{R} - \Delta SWC - DP - E_{S}$$

(1)

where, T_wb corresponds to daily transpiration obtained by the water balance, P is precipitation, I_R is the amount of water applied through irrigation, ΔSWC is the difference in soil water content between two consecutive days, DP is deep percolation and E_S corresponds to evaporation. However, the water balance was calculated during days without P and no I_R applied, therefore P, DP and I_R were considered null. During sap flow validation, the soil was covered with plastic sheeting to avoid evaporation fluxes (E_S = 0). Soil water content (SWC) was measured every 20 cm down to 180 cm depth using a neutron probe (Campbell Pacific Nuclear Scientific, Model 503). The neutron probe measurements were calibrated based on the volumetric moisture content (cm³ of water/cm³ of soil) of soil samples taken at the time of tube installation. Tube design consisted of six tubes in each tree installed in one quarter of the planting area. Two groups of three tubes (6 tubes) were installed in parallel. The distribution of each group was: below the emitter, at a quarter of the distance between rows and at half the distance between rows. Calibration coefficients estimated with T_wb and T_sf measurements were assumed to be constant throughout the season (Espadafor et al. 2015). The calibration coefficient varied from tree to tree and ranged from 0.56 to 1.62. The calibrated T_sf was used to calculate total hourly transpiration at the time of image acquisition.

Stem water potential, leaf area index and fractional intercepted photosynthetically active radiation

In those trees with sap flow data, the midday stem water potential (Ψ_stem) was measured continuously every two weeks following the protocol established by McCutchan and Shackel (1992). The Ψ_stem was calculated as an average of three measurements taken from each tree on ten dates during the growing season. Shaded leaves were selected and kept in a plastic bag covered by aluminum foil for 1 h before the measurement to equilibrate the water potential between leaf, stem and branches. All measurements were acquired in less than 1 h with a pressure chamber (Plant Water Status Console, Model 3500; Soil Moisture Equipment Corp., Santa Barbara, CA).

The LAI was measured using both direct and indirect methods. The direct measurement corresponded to a destructive sampling in eight trees, one for each production system in the fully irrigated treatment and on two different dates corresponding to March 25th and June 8th of 2021. The total amount of leaf extracted was weighed and the leaf area of 200 g of sample was measured using the LI-3100C Area Meter instrument (LI-COR Inc., Lincoln, NE, USA). The total LAI was extrapolated based on a ratio of sample weight to leaf area.

The LAI-2200 Plant Canopy Analyzer (PCA) (LI-COR Inc., Lincoln, NE, USA) was used to obtain LAI indirectly in trees with sap flow measurements and for each flight campaign at midday. The radiation from the sky incident above the tree was measured in an open space using five sensor’s rings. Subsequently, four measurements were made under the tree in the N, S, E and W directions. Measurements were performed by covering the sensor lens with the diffuse cap. The data were imported into the FV2200 v. 2.1.1 software to calculate the LAI using the vertical profile of the tree crown. The vertical profile was obtained using canopy height images of a digital elevation model (DEM) generated from point clouds of multispectral images.

In addition, the diurnal evolution of fractional intercepted photosynthetically active radiation (fIPAR) was obtained for all measured trees in each production system following the methodology described by Casadesús et al. (2011). The fIPAR was obtained by installing in each tree 125 cm long sensor bars containing photodiodes every 10 cm (VTB8440BH, PerkinElmer Optoelectronics, Vaudreuil, Canada). Each photodiode was installed in an aluminum channel with an outer section of 25 mm and covered by a 4.6 mm thick PTFE sheath (TECAFLON PTFE, Ensinger Ltd., Llantrisant, UK) that acted as a light diffuser. The photodiodes have a spectral response over the photosynthetically active radiation (PAR) region between wavelengths 330 and 720 nm (with a peak at 580 nm). The PAR below the canopy (PAR_below) was measured using bars installed on the ground of the planting system. The sensing bars were placed parallel to the crop rows and positioned to cover all the tree spacing. Due to the different spacing between production systems, 10 sensing bars were installed in the open vase system, while 5 were used for the central axis and hedgerow systems. For its part, PAR above the canopy (PAR_above) was obtained from two sensing bars installed outside the field with direct incident PAR. The PAR data were recorded every 15 min from March to October of 2021. Then, the average hourly fIPAR at the time of image acquisition was calculated.

Image acquisition campaign

The image acquisition campaign consisted of five flights conducted on March 24, May 19, June 3 and 29, and July 29 of 2021. All flights were carried out with an UAV Dronehexa XL (DRONETOOLS, Seville, Spain) equipped with a multispectral and thermal camera. The Micasense RedEdge-MX (Micasense, Northlake Way, Seattle, USA) has five spectral bands located at the wavelengths 475 ± 20 nm, 560 ± 20 nm, 668 ± 10 nm, 717 ± 10 nm, and 840 ± 40 nm. The thermal camera used was the FLIR SC655 (FLIR Systems, Wilsonville, OR, United States) which has a spectral response in the range of 7.5–13 µm and an image resolution of 640 × 480 pixels. All flights were conducted at solar time (14:00h local time) under clear sky conditions and at 50 m above ground level to obtain images with a spatial resolution of 0.03 m and 0.06 m for the multispectral and thermal images, respectively.

Images were radiometrically, atmospherically and geometrically corrected. The radiometric calibration for the multispectral images was applied using an external incident light sensor which measured the irradiance levels of light at the same bands as the Micasense multispectral sensors. In addition, in situ spectral measurements for ground calibration targets were performed using a Jaz spectrometer (Ocean Optics, Inc., Dunedin, FL, United States) for radiometric calibration. The Jaz has a wavelength response from 200 to 1,100 nm and an optical resolution of 0.3–10.0 nm. During spectral collection, spectrometer calibration measurements were taken with a reference panel (white color SpectralonTM) and dark current before and after taking readings from radiometric calibration targets. The radiometric calibration of the thermal sensor was assessed in the laboratory using a blackbody (model P80P, Land Instruments, Dronfield, United Kingdom). In addition, in situ measurements were conducted in the field concomitant to image acquisition in different ground calibration targets. In situ temperature measurements were conducted with an SI-111-SS apogee infrared radiometer connected to an Apogee AT-100 microCache Bluetooth micrologger (Apogee instruments Inc, Logan, UT, USA). The mosaicking process of obtaining thermal and multispectral images and the generation of the DEM from point clouds of multispectral images was performed with the Agisoft Metashape Professional software (Agisoft LLC., St. Petersburg, Russia). After the mosaicking process, geometric and radiometric corrections were carried out using QGIS 3.4 (QGIS 3.4.15).

Temperatures and biophysical traits using UAVs

To retrieve T_c and T_s separately and the biophysical traits, a canopy layer was obtained. The canopy layer was created based on a contextual classification using the DEM and the soil-adjusted vegetation index (SAVI). SAVI minimizes the influence of soil brightness in the red and near-infrared wavelengths, thereby enhancing the contrast between vegetation and the soil surface (Qi et al. 1994). Pixels with a DEM greater than 1.5 m and SAVI greater than 0.2 were classified as canopy, while pixels that did not meet the conditions were classified as soil. The canopy layer was used to extract the T_c and T_s from thermal images, while T_rad corresponds to the total average temperature in each scene. The fractional canopy cover (f_c), canopy height (h_c) and canopy width (w_c) were obtained using the canopy layer and DEM (Fig. 2). In addition, the canopy volume (v_c), normalized difference vegetation index (NDVI), normalized difference water index (NDWI) and modified triangular vegetation index (MTVI2) were calculated as additional inputs of a machine learning approach for estimating LAI in all trees. The extraction of biophysical traits and temperature from high-resolution images was conducted using the Python programming language (Python Software Foundation. Python Language Reference, version 3.10. Available at http://www.python.org).

Fig. 2

Flowchart of the procedures used for processing the multispectral and thermal images and the digital elevation model (DEM) to obtain the different biophysical variables of the vegetation and some of the inputs needed in the different two-source energy balance (TSEB) modelling approaches

According to Gao et al. (2022), who compared machine learning algorithms to estimate LAI, the random forest technique performed slightly better than the other algorithms. Therefore, the random forest algorithm (scikit-learn Python library) was trained in this study to estimate LAI using the contextual, spectral, structural information as input data (production system, v_c, f_c, w_c, h_c and the canopy mean and canopy standard deviation of the vegetation indexes NDVI, NDWI and MTVI2) and the measured LAI as calibration data. The model was calibrated using a random sample of 80% of the total collected data. Figure 2 shows the methodological scheme used to obtain the main parameters needed for the TSEB model.

TSEB model description

Transpiration was retrieved from the TSEB model originally formulated by (Norman et al. 1995) and further improved by Kustas and Anderson (2009). The TSEB model is based on an energy balance approach that assumes that the energy available at the surface is distributed mainly between sensible heat flux (H), latent heat flux (LE) and soil heat flux (G). Therefore, LE (W m⁻²) can be calculated as a residual of the surface energy equation (Eq. 1):

$$LE \approx R_{n} - H - G$$

(1.1)

$$LE_{s} \approx R_{n, s} - H_{s} - G$$

(1.2)

$$LE_{c} \approx R_{n, c} - H_{c}$$

(1.3)

where R_n is net radiation and the subscripts _C and _S correspond to the canopy and soil covers, respectively. The R_n,s and R_n,c were calculated using the canopy radiative transfer model of Campbell and Norman (1998). Radiation models depend on accurate estimation of shortwave transmittance (τ_C). The radiative transfer model of Campbell and Norman (1998) has been widely used to estimate τ_C in energy balance models. However, woody crops are usually arranged in rows, so τ_C models have to add the effect of radiation transmitted to the surface through the inter-row space to the radiation transmitted through the canopy (Kustas and Norman 1997, 1999a, b; Anderson et al. 2005). Due to their accuracy and simplicity, two models were selected to estimate τ_C: (1) a basic clumping index meant for heterogeneous randomly placed clumped canopies combined with the Campbell and Norman transfer model (C&N–H); and (2) a rectangular hedgerow clumping index combined with the Campbell and Norman transfer model (C&N–R). The dependency of the radiative transfer τ_C on wavelength is accounted for, in part, by a leaf absorption parameter for the PAR and near infrared spectra, which may be species dependent (e.g., Gausman and Allen 1973). The estimated PAR transmittance (τ_C,PAR) was evaluated using fIPAR obtained from the sensing bars, where fIPAR corresponds to 1-τ_C,PAR. More detailed information about the two transfer models can be found in Parry et al. (2019), who evaluated and compared the radiation partitioning of both models in a vineyard.

G is usually calculated as a constant fraction of R_n,s of around 0.35 at solar noon. H is partitioned into soil (H_S) and canopy (H_C) fluxes, with their corresponding resistances acting in series—as an analogy of Ohm’s law for electrical transport (Eqs. 2.1, 2.2, 2.3).

$$H_{s} = \rho C_{p} \frac{{T_{s} - T_{ac} }}{{R_{s} }}$$

(2.1)

$$H_{c} = \rho C_{p} \frac{{T_{c} - T_{ac} }}{{R_{X} }}$$

(2.2)

$$H_{s} + H_{c} = \rho C_{p} \frac{{T_{ac} - T_{a} }}{{R_{A} }}$$

(2.3)

where $\rho$ is the air density, ${C}_{p}$ is the specific heat of air, T_s is the soil temperature, T_c is the canopy radiometric temperature, T_a is the air temperature, T_ac is the temperature in the canopy air space, R_s is the resistance to heat flow in the boundary layer immediately above the soil surface, R_x is the total boundary layer resistance of the complete canopy leaves, and R_a is the aerodynamic resistance to turbulent heat transport between the air canopy layer and the overlying air layer. The resistances were estimated following Kustas and Norman (1997, 1999a, b) and Norman et al. (1995).

Energy fluxes can be directly estimated by contextual estimation of canopy and soil temperature (TSEB-2T) using high resolution temperature images, as shown in Fig. 2a. However, there are cases where it is difficult to obtain these temperatures separately. For this reason, the Priestley-Taylor iterative retrieval model (TSEB-PT) was developed by Norman et al. (1995) as a solution to estimate T_c and T_s by a single directional radiometric temperature observation T_rad using Eq. 3.

$$T_{rad} \left( \theta \right) \approx [f_{c} \left( \theta \right)T_{c}^{4} + \left( { 1 - f_{c} \left( \theta \right) } )T_{s}^{4} \right]^{0.25}$$

(3)

where ${f}_{c}\left(\theta \right)$ corresponds to the fraction of vegetation observed at a zenith angle $\theta$. An iterative process was established to estimate H_S, H_C, T_c and T_s, defining an initial guess of the potential canopy transpiration.

$$LE_{c} = \alpha_{PT} f_{c} \frac{{\Delta }}{{{\Delta } + \gamma_{p} }}R_{n, c}$$

(4)

where ${\alpha }_{PT}$ is the Priestley-Taylor coefficient initially set to 1.26, f_C is the fraction of LAI that is green and thus capable of transpiring, $\Delta$ is the slope of saturation vapor pressure vs. temperature, and ${\gamma }_{p}$ is the psychometric constant. Using Eqs. 2, 3 and 4, LE_C is iteratively recalculated by reducing the ${\alpha }_{PT}$ value until both LE_S and LE_C become zero or positive. The calculated LE_C was converted to hourly transpiration in mm using: 1000 $\times$ 3600 $\times$ LE_C/($\rho$ _w λ) where 1000 converts m to mm, 3600 converts seconds to hours, $\rho$ _w is the density of water (assumed to be 1000 kg m⁻³) and λ is the latent heat of vaporization (J kg⁻¹): λ = 1e⁶ $\times$ (2.501–0.002361 T_a).

This study compares the TSEB-2T and TSEB-PT approaches to estimate almond transpiration in combination with the Campbell and Norman transfer model with a basic (C&N–H) and with a rectangular hedgerow clumping index (C&N–R) to estimate τ_C. Therefore, four models will be assessed: TSEB-PT_H, TSEB-PT_R, TSEB-2T_H and TSEB-2T_R. The subscript _H represents the C&N–H shortwave transmittance model and the subscript _R represents the C&N-R shortwave transmittance model. TSEB models developed by Nieto et al. (2019), which are available in the Python programming language, were used (https://github.com/hectornieto/pyTSEB). Evaluation of the TSEB transpiration models is performed through a comparison with hourly T_sf measurements at the time of image acquisition for each of the five flights performed during the 2021 season.

Model evaluation

Model agreement was evaluated using the following metrics: correlation coefficient (R², Eq. 5), bias (Eq. 6), Root Mean Square Error (RMSE, Eq. 7) and Mean Absolute Error (MAE, Eq. 8).

$$R^{2} = \,\frac{{\sum\nolimits_{i = 1}^{n} {\left( {E_{i} - \overline{E} } \right)\,\left( {O_{i} - \overline{O} } \right)} }}{{\left[ {\sum\nolimits_{I = 1}^{N} {\,\left( {O_{i} - \overline{O} } \right)^{2} } } \right]^{0.5} \left[ {\sum\nolimits_{I = 1}^{N} {\left( {E_{i} - \overline{E} } \right)^{2} \,} } \right]^{0.5} }}$$

(5)

$$ bias= \frac{{\mathop \sum \nolimits_{i = 1}^{n} \left( {E_{i} - O_{i} } \right)}}{n}$$

(6)

$$RMSE = \sqrt {\frac{{\mathop \sum \nolimits_{i = 1}^{n} \left( {E_{i} - O_{i} } \right)^{2} }}{n}}$$

(7)

$$MAE = \frac{{\mathop \sum \nolimits_{i = 1}^{n} \left| {O_{i} - E_{i} } \right|}}{n}$$

(8)

$$NASH = 1.0 - \frac{{\mathop \sum \nolimits_{i = 1}^{n} \left( {O_{i} - E_{i} } \right)^{2} }}{{\mathop \sum \nolimits_{i = 1}^{n} \left( {O_{i} - \overline{O}} \right)^{2} }}$$

where $n$ is the number of observations, ${O}_{i}$ is the measured value, ${E}_{i}$ is the estimated value, $\overline{O }$ is the average of measured values and $\overline{E }$ is the average of estimated values. The evaluation of the model was performed using the Python programming language.

Results

Biophysical traits

The actual LAI of defoliated trees showed R² and RMSE values of 0.89 and 0.57 m² m⁻², respectively, with respect to estimated LAI using the LAI-2200 Plant Canopy Analyzer (LAI₂₂₀₀) (Fig. 3a). Trees with a higher LAI corresponded with those defoliated on June 8, while those with a lower LAI corresponded with those defoliated on May 25. Eleven variables were required for the random forest model. These correspond to vegetation indices obtained from multispectral images (mean MTVI, std MTVI, std NDWI, mean NDWI, std SAVI, mean NDVI), structural parameters of trees (f_c, v_c, w_c, h_c) and production system as contextual information. When the random forest LAI (LAI_modeled) was regressed against LAI₂₂₀₀, it showed R² and RMSE values of 0.94 and 0.30 m² m⁻², respectively (Fig. 3b). It can be observed that there was high variability in LAI values, which ranged from 1.4 to 4.5 m² m⁻². LAI estimations presented significant differences in both production system and irrigation treatment, but not in their interaction (Table 1). The hedgerow system had a mean LAI of 3.42 m² m⁻² during the season and was significantly higher than the open vase and central axis systems (3.19 and 3.03 m² m⁻², respectively) (Table 2). The fully irrigated treatment showed significantly higher values of LAI, averaging 3.33 m² m⁻², while the severe stress treatment presented the lowest, averaging 3.10 m² m⁻².

Fig. 3

Calculated vs. measured leaf area index (LAI) for LAI₂₂₀₀ vs. LAI defoliated (a) and LAI estimated by random forest using inputs derived from UAV images (LAI_modeled) vs. LAI₂₂₀₀ (b)

Table 1 Analysis of variance (three-way ANOVA) testing the influence of date, production system (PS) and irrigation treatment (TRT) on biophysical traits during the flight campaign

Table 2 Comparison of the main biophysical traits measured during the flight campaign

Other evaluated biophysical traits showed significant differences in both production system and irrigation treatment, with the partial exception of f_c for the latter (Table 1). The fIPAR measured at the time of image acquisition did not show significant differences in either production system or irrigation treatment. However, when considering the interaction between production system and irrigation treatment, the trees measured in the hedgerow system with the fully irrigated treatment presented significantly lower values of fIPAR. Table 2 shows the mean values for all the biophysical traits in each production system (mean row), the interaction between production system and irrigation treatment, and also Tukey’s analysis. In general, the open vase systems had larger trees with higher h_c, and f_c, while the hedgerow system had higher values in h_c/w_c than the open vase and central axis systems. The trees from the severe stress irrigation treatment were slightly smaller throughout the season (h_c of 4.20 m, f_c of 0.52) than those from the mild stress (h_c of 4.67 m, f_c of 0.57) and fully irrigated treatments (h_c of 4.56, f_c of 0.54).

Water applied and physiological measurements

The total amount of irrigation water applied (IR) in each production system and each irrigation treatment is shown in Table 2. Differences in the total amount of water applied were statistically significant between irrigation treatments (p < 0.0001) but not between production systems. Stem water potential (Ψ_stem) showed significant differences between production system, irrigation treatment and their interaction (p < 0.0001) (Table 1). This is mostly attributable to the open vase system which, for the latest growing stages, showed significantly higher values in the mild stress and severe stress treatments (Table 2).

The seasonal pattern of Ψ_stem by production system and irrigation treatment is shown in Fig. 4. At the beginning of the growing season, all trees had a Ψ_stem around − 0.6 MPa. No significant differences were observed between irrigation treatments in this initial growing stage. After May 5, treatments began to differentiate and Ψ_stem decreased until June 15 when the values stabilized. During this period, the fully irrigated and mild stress treatments showed Ψ_stem values above − 1.0 MPa, without significant differences between them. However, the Ψ_stem value of the severe stress treatment declined progressively until reaching − 1.6 MPa. From June 15 to harvest, the severe stress treatment had significantly lower values for all the production systems, ranging between − 1.4 and − 2.1 MPa. During this period, the Ψ_stem values of the open vase system were slightly higher in comparison to those for the central axis and hedgerow systems. The open vase system had values ranging from − 1.0 to − 1.7 MPa, while the central axis and hedgerow values ranged between − 1.1 and − 2.1 MPa. On the other hand, the Ψ_stem exhibited significant differences between fully irrigated and severe stress treatments for all dates, but some similarity between the fully irrigated and mild stress treatments, except in the central axis and hedgerow systems on July 29.

Fig. 4

Seasonal trends of midday stem water potential (Ψ_stem) in each production system and irrigation treatment. Dashed lines indicate the days on which flights were performed (color figure online)

Transpiration measured with sap flow sensors (T_sf) also showed significant differences in production system and irrigation treatment, but not in their interaction (Table 1). The seasonal pattern of daily transpiration by production system and irrigation treatment is shown in Fig. 5. This clearly shows higher transpiration values in the fully irrigated treatment and lower values in the severe stress treatment for all the production systems. Transpiration began to differentiate between treatments around ten days after May 5 when irrigation was applied. Thus, no significant differences in T_sf were observed for the first flight conducted on March 24. On that date, T_sf ranged from 0.86 to 1.78 mm day⁻¹. Transpiration in the fully irrigated treatment followed a similar pattern to ETo, with the maximum values observed during July coinciding with dates with higher ETo and high canopy development. These maximum values reached 7.8 mm day⁻¹ for some specific dates. The maximum values of transpiration for the severe stress treatment were observed during May and June, depending on the production system. During August and September, T_sf started progressively declining in all the treatments with ETo and leaf senescence. In the open vase and hedgerow systems, transpiration in the mild stress and fully irrigated treatments followed a similar pattern of values throughout the season. However, for the central axis system, the values for the mild stress treatment were closer to those for the severe stress treatment.

Fig. 5

Seasonal trends of daily sap flow transpiration (T_sf) measurements during 2021 for each production system and irrigation treatment. Dashed lines correspond to reference evapotranspiration (ETo) (color figure online)

Table 2 shows the average hourly transpiration measured at the time of image acquisition. The open vase system had significantly higher transpiration rates in all treatments. Overall, cumulative transpiration in the fully irrigated treatment was 761, 625 and 487 mm for the open vase, central axis and hedgerow systems, respectively. The mild stress treatment accumulated 701, 545, and 520 mm of transpiration for the open vase, central axis and hedgerow systems, respectively. Finally, the cumulative transpiration of the severe stress treatment was 594, 495, 384 mm for the open vase, central axis and hedgerow systems, respectively. In the case of the hedgerow system, trees of the mild stress irrigation treatment transpired 33.9 mm more than those of the fully irrigated treatments during the season. This is explained by its larger canopy size and similar Ψ_stem values until July 29, just before ETo began to decrease. In Fig. 5, symbols in bold correspond to T_sf values recorded for the flight campaign dates. It can be observed that, except for the first flight (March 24), all flight dates showed notable differences between irrigation treatments in each production system.

Comparison of the Priestley-Taylor and contextual method for the retrieval of T _c and T _s

The comparison of the T_c obtained with the contextual methodology (T_c2T) and the two Priestley-Taylor methodologies, applying the shortwave transmittance models C&N–H (T_cPT_H) and C&N-R (T_cPT_R), showed an overall RMSE of 1.46 °C and 1.24 °C, respectively (Fig. 6). Thus, the adoption of the C&N-R model resulted in an estimate of T_c closest to the contextual method. However, differences in T_c between the contextual and the T_cPT_H model, and between the contextual and the T_cPT_R model varied significantly between dates, production system, irrigation treatment and the interaction of the latter two (p < 0.0001). The T_cPT_H and T_cPT_R models had significantly higher errors in the hedgerow production system (RMSE of 2.03 and 1.45 ºC for T_cPT_H and T_cPT_R, respectively). The T_cPT_H and T_cPT_R showed a similar error in the open vase and central axis systems (RMSE values of 0.95 and 1.04 °C in the open vase system for T_cPT_H and T_cPT_R, respectively, and 1.21 and 1.21 °C in the central axis system for T_cPT_H and T_cPT_R, respectively) (data not shown). In addition, the T_cPT_H and T_cPT_R approaches had significantly higher (p < 0.0001) errors in the fully irrigated treatment (RMSE of 1.68 and 1.37 °C for T_cPT_H and T_cPT_R, respectively) in comparison to the other treatments. This was especially evident in the hedgerow system, which had a higher error in both models.

Fig. 6

Comparison of canopy temperature (T_c) obtained with the contextual (T_c2T) approach and T_c obtained with the Priestley-Taylor approach using the C&N–H (T_cPT_H) (a) and the C&N-R (T_cPT_R) (b) transmittance models. The color indicates each irrigation treatment (TRT), with blue, orange and red representing the fully irrigated (FI), mild stress (MS) and severe stress (SS) treatments, respectively. The point shape shows each production system (PS) with OV, CA and HGR corresponding to open vase, central axis and hedgerow, respectively (color figure online)

In addition, T_cPT_H and T_cPT_R showed significantly lower errors in the mild stress treatment. Figure 7 shows a comparison of the T_c–T_a mean retrieved with the PT_H, PT_R and contextual approaches. Overall, the PT_H model estimated significantly higher T_c–T_a values than the PT_R and contextual approaches. Moreover, significant differences between models were observed in the hedgerow system, but not in the open vase and central axis systems. The PT_H model overestimated the T_c–T_a in the hedgerow system with a mean 3.53 ºC vs. 2.65 ºC and 2.44 ºC in the PT_R and contextual approaches, respectively.

Fig. 7

Differences between irrigation treatments (TRT) and productivity systems (PS) in the average T_c-T_a retrieved with the Priestley-Taylor approach using the C&N–H (T_cPT_H) and the C&N-R (T_cPT_R) transmittance models and the contextual (T_c2T) approach. Mean columns correspond to the average value of all irrigation treatments. Different letters indicate significant differences at p < 0.05 using Tukey’s test (color figure online)

Assessment of the radiative transfer models for the retrieval of fIPAR

Modeled fIPAR with the C&N–H and C&N-R approaches was regressed with measured fIPAR (Fig. 8a, b). The C&N–R model showed a better agreement than the C&N–H model, with an RMSE of 0.09 and 0.12, respectively. Both the C&N–H and C&N–R models overestimated fIPAR with a bias of 0.09 and 0.03, respectively, which represented an error of 17.5% and 6.1% with respect to measured fIPAR. However, C&N–H estimated significantly higher values than C&N-R and measured fIPAR (0.59 vs. 0.53 and 0.50 for C&N–H vs. C&N–R and measured fIPAR).

Fig. 8

Relationship between estimated and measured hourly fractional intercepted photosynthetically active radiation (fIPAR) for the canopy radiative transfer model C&N–H (a) and C&N–R (b)

The accuracy of the fIPAR estimates depended on the production system (Fig. 9). In the hedgerow production system, there was better agreement with the C&N-R (RMSE of 0.074) than with the C&N–H (RMSE of 0.156). The latter estimated significantly higher values than both C&N-R and measured fIPAR. In the central axis system, C&N-R showed a lower error than C&N–H (RMSE of 0.69 vs. 0.091 in C&N–H vs. C&N–R). In contrast, the fIPAR with the model C&N-R in the open vase production system indicated a higher error in comparison to the C&N–H (RMSE of 0.095 vs. 0.125 for C&N–H vs. C&N-R). Open vase trees present a more isolated structure than those in central axis and hedgerow systems, which explains the better fit with the C&N–H model. Despite this, significant differences between models were not observed in either the open vase or central axis systems.

Fig. 9

Comparison between the mean fIPAR values obtained in the two radiative transfer models with the measured values for each production system (color figure online)

Validation of transpiration with different TSEB methods

Table 3 shows the fluxes retrieved with the different TSEB modeling approaches for each production system and irrigation treatment. Only the central axis and hedgerow systems showed significant differences between models in R_n,c, H_C and LE_C. Models that used the C&N–H approach (TSEB-PT_H and TSEB-2T_H) estimated R_n,c values significantly higher (p < 0.0001) than the others. This occurred for all irrigation treatments. For H_C, TSEB-PT_H indicated significantly higher values (p < 0.0001) than the other models, which did not show significant differences between them. Similarly, the TSEB-2T_H model showed significantly higher values of LE_C (p < 0.0001) than the others, between which no significant differences were observed. Estimated LE_S ranged from 0 to 17.5% of total LE and was significantly different between models, production systems and irrigation treatments. The TSEB-PT_R and TSEB-2T_R estimated significantly higher LE_S than TSEB-PT_H and TSEB-2T_H. With respect to the production systems, the hedgerow system presented higher estimated LE_S, while the central axis system recorded significantly lower LE_S values. For its part, significantly lower LE_S values were estimated in the severe stress treatment.

Table 3 Estimated energy partitioning into canopy net radiation (R_n), canopy latent heat flux (LE_C), canopy sensible heat flux (H_C), soil latent heat flux (LE_S), soil sensible heat flux (H_C) and soil heat flux (G), retrieved with the different TSEB modelling approaches for each production system and irrigation treatment

In all production systems, the TSEB-2T_R model had the best performance in terms of almond crop transpiration, with an RMSE of 0.13 mm h⁻¹ and a bias of 0.08 mm h⁻¹ when regressed against T_sf (Fig. 10 and Table 4). In contrast, the TSEB-PT_H, TSEB-PT_R and TSEB-2T_H models showed RMSE values of 0.22, 0.17 and 0.21 mm h⁻¹, respectively. We observed that use of the C&N-R transmittance model improved the estimates of transpiration in both the TSEB-PT and TSEB-2T models when compared to C&N–H. These improvements were reflected in the RMSE (0.22 vs 0.17 mm h⁻¹ for TSEB-PT_H vs. TSEB-PT_R and 0.21 vs 0.13 mm h⁻¹ for TSEB-2T_H vs TSEB-2T_R). Adoption of the C&N–R transmittance model also contributed to a significant reduction in R², bias and the mean absolute error (MAE) in both TSEB-PT and TSEB-2T. Moreover, the hedgerow system presented significantly higher RMSE in all the TSEB models, except in TSEB-PT_R (Table 4). Regarding irrigation treatments, only the TSEB-PT_H model did not show significant differences in the RMSE between irrigation treatments. The TSEB-PT_R showed a higher RMSE in the severe stress and mild stress treatments than in the fully irrigated treatment. For their part, the two TSEB-2T models showed a higher RMSE in the fully irrigated treatment than in the mild stress and severe stress treatments.

Fig. 10

Regressions between measured and estimated hourly transpiration with the modelling approaches TSEB-PT_H (a), TSEB-PT_R (b), TSEB-2T_H (c), and TSEB-2T_R (d). The color indicates each irrigation treatment (TRT), with blue, orange and red representing the fully irrigated (FI), mild stress (MS) and severe stress (SS) treatments, respectively. The point shape shows each production system (PS), with OV, CA and HGR corresponding to open vase, central axis and hedgerow, respectively (color figure online)

Table 4 Root mean squared error (RMSE) and significant differences between hourly estimated transpiration rates, using the different TSEB models vs. transpiration measured with sap flow sensors, for each production system (PS) and irrigation treatment (TRT)

Discussion

Many studies have reported different values for the almond crop coefficient (K_C), with these mostly depending on environmental conditions, water management, fractional canopy cover, fIPAR and albedo (Girona 2005, Espadafor et al. 2015; García-Tejero et al. 2015). In our results, the T_sf had significant differences between production systems, which may be affected by h_c, f_c, h_c/w_c, and LAI. Considering K_T as a ratio between transpiration and reference evapotranspiration (T/ETo) in fully irrigated conditions, the maximum K_T values obtained in this study were 1.2, 1.13 and 0.94 for the open vase, central axis and hedgerow production systems, respectively. These values are slightly higher than those reported in other recent studies. For instance, Bellvert et al. (2018) estimated K_T values of 1.0 in almond trees planted at a spacing distance of 5.5 × 7.3 m and with an LAI of around 3 m m⁻². Espadafor et al. (2015) also estimated a maximum K_T of 1.02 in almonds with an f_c of around 50% based on a ratio with daily fIPAR. In another study, Espadafor et al. (2017) reported a K_T of around 1.1 in almonds with an f_c of around 40%. In our study, however, the K_T values corresponded to almonds with an f_c higher than 55% and LAI values between 3.5 and 4.5 m m⁻². Therefore, the higher size of trees in the open vase and central axis systems may explain why the K_T values obtained in this study were higher than others reported in the literature. The fact that fully irrigated trees in the hedgerow system had the lowest K_T can be explained by their smaller size (Table 2). Differences between canopy architectures have a direct effect on the amount of light intercepted by the canopy and consequently on the transpiration rates. We observed a variability in transpiration rates of 0.11–0.83 mm h⁻¹, from the beginning of the season in March to the maximum transpiration values in July, with significant differences among production systems and irrigation treatments.

Our results showed an RMSE error of 0.30 m² m⁻² in UAV-modeled LAI, which corresponds to a relative RMSE of around 11%. These results are slightly better than those reported by Gao et al. (2022) who obtained an error of 0.31 m² m⁻² (25.3% of relative RMSE) using random forest algorithms. In contrast to the latter, our study also incorporated in the random forest model canopy structure and contextual information of trees under different production systems. Bellvert et al. (2021) reported, for almond trees, an RMSE of 0.24 m² m⁻² using a multiple regression model with parameters of canopy structure (crown area and canopy volume) estimated from multispectral images. However, the lower error recorded by Bellvert et al. (2021) can be attributed to the fact that the trees had a lower LAI (between 0.5 and 2 m² m⁻²) compared to the values obtained in this study (between 1.4 and 4.5 m² m⁻²). Gao et al. (2022) demonstrated that use of a hybrid machine learning approach with random forest and relevance vector machine algorithms provides better results than random forest by itself. This type of approach could be easily applied in future studies using our datasets to assess the degree of improvement in LAI estimates.

It was observed that regression between transpiration estimations with the different TSEB approaches and sap flow measurements had overall RMSE values of 0.22, 0.18, 0.21 and 0.13 mm h⁻¹ for the TSEB-PT_H, TSEB-PT_R, TSEB-2T_H and TSEB-2T_R models, respectively. These values correspond to relative RMSE values of 61%, 49%, 56% and 36%, respectively. Numerous studies have evaluated different approaches of the TSEB model in woody crops, including grapevine, olive and almond (Nieto et al. 2019; Aguirre-García et al. 2021; Guzinski et al. 2021; Kool et al. 2021; Nassar et al. 2021; Jofre-Čekalović et al. 2022; Peddinti and Kisekka 2022). However, to our knowledge, no study has assessed ET partitioning and validated the transpiration component with sap flows in trees with different production systems and irrigation treatments.

Those studies that validated the transpiration component were based on different modelling approaches with eddy covariance flux tower (Kool et al. 2021; Zhang et al. 2022; Gao et al. 2023; Knipper et al. 2023) data and sap flow measurement (Peng et al. 2023). Among them, Kool et al. (2021) obtained the highest accuracies with the TSEB-PT approach, indicating an error of 35 W m⁻² (29% relative error) in grapevines. In contrast, Gao et al. (2023) reported an error of 70 W m⁻² and an R² of 0.54 using the TSEB-2T model, also in grapevines. In almonds, Knipper et al. (2023) obtained an error of 0.82 mm day⁻¹ with the TSEB-PT using satellite data at 30 m resolution. One of the main challenges faced in these works was accurately measuring actual transpiration. To address this issue, Gao et al. (2023) compared three eddy covariance flux tower methods for ET partitioning and found that the Modified Relaxed Eddy Accumulation (MREA) and Conditional Eddy Covariance (CEC) produced the most consistent results. However, these eddy covariance models rely on various assumptions and their outcomes are estimates, rather than precise transpiration values. Furthermore, they are unable to differentiate between crop and inter-row transpiration. To address this uncertainty, the synergistic application of sap flow and water balance methods has proven effective for monitoring transpiration in orchards (López-López et al. 2018a; Peng et al. 2023). For instance, Peng et al. (2023) successfully established a reliable relationship between transpiration estimates using TSEB-PT and actual transpiration measured through sap flow calibrated with water balance methods in tomato crops. Nevertheless, measuring transpiration in complex woody crops, such as almond, remains a challenge. It is therefore crucial to acknowledge that the accurate measurement of actual transpiration in almond cultivation continues to present a challenge.

In our study, the percentage of estimated evaporation with respect to the total estimated ET represented an average of 6.68% and ranged from 0 to 17%. Due to the low estimated evaporation, it is possible to perform a relative comparison of the transpiration evaluated in this study with previous evaluations of ET or LE in almond crops (Xue et al. 2020; Jofre-Čekalović et al. 2022; Peddinti and Kisekka 2022). In their study, Xue et al. (2020) reported RMSE values of 0.93, 1.19 and 1.53 mm day⁻¹ for the Surface Energy Balance Algorithm for Land (SEBAL), Mapping Evapotranspiration at High Resolution with Internalized Calibration (METRIC) and Surface Energy Balance System algorithm (SEBS) models, respectively. For instance, Peddinti and Kisekka (2022) reported an relative RMSE of 10% (0.64 mm day⁻¹ of RMSE) using the TSEB-PT approach. Jofre-Čekalović et al. (2022) also observed an average relative RMSE of 30% (equivalent to 87 W m⁻² of RMSE) using the sharpened temperature (TSEB_S2+S3) and TSEB-PT with almonds that exhibited a high degree of water status heterogeneity. In our study, we obtained slightly higher errors which, in part, may be attributed to the fact that our validations were conducted with in situ sap flow sensors; this gave us actual measurements rather than estimations. Additionally, our evaluation included trees with a higher degree of variability in their canopy architecture and water status than in other studies. We also evaluated a greater number of trees. Among the different shortwave transmittance models evaluated, the lowest error was obtained with the C&N-R. This emphasizes the need to use radiative transfer models adapted to the canopy shapes of woody row crops. When comparing the R_n,c values obtained with the two transmittance models, it can be observed that the C&N-R model estimated a radiation absorbed by the canopy which was 3%, 10% and 18% lower than the C&N–H model for the open vase, central axis and hedgerow production systems, respectively. Similarly, Parry et al. (2019) reported an underestimation of transmittance (R_n,c is inversely proportional to transmittance) using the C&N–H model in comparison to the C&N–R model. The same study showed a higher sensitivity of the C&N–H model to LAI. In our study, despite the high LAI values in the hedgerow production system, the amount of radiation transmitted in the inter-row space was of greater importance due to the lower f_c and higher h_c/w_c ratio. Therefore, we hypothesize that the highest values of LAI observed in the hedgerow system probably affected the transmittance estimations with the C&N–H model. Along with the R_n,c reduction, the C&N-R model showed a better agreement in T_c estimates than the C&N–H when using the PT approach (1.46 vs. 1.24 °C for T_cPT_H vs. T_cPT_R). This fact explains the better accuracy of the TSEB-PT_R modelling approach to estimate transpiration. These differences in the TSEB-PT models are in concordance with Kustas and Norman (1997, 1999a, b) and Anderson (2005), who demonstrated a reduction in overestimation of the sensible heat fluxes when TSEB-PT was used considering the clumping effect in the transmittance models.

The TSEB-2T and C&N-R transmittance models synergistically (TSEB-2T_R) obtained the most robust estimations of canopy transpiration. In grapevine, Nieto et al. (2019) concluded that using a contextual algorithm to derive soil and canopy temperatures separately yielded a closer agreement with flux tower measurements in comparison to the Priestley-Taylor model. Gao et al. (2023) also reported better agreement in ET partitioning using a temperature-separation method named quantile temperature separation (QTS). In contrast, Kool et al. (2021) obtained much lower accuracies using the TSEB-2T model. They attributed this to the loss of flexibility within the model or a bias in measured temperatures. In any case, improving on techniques based on the combination of photogrammetry, multispectral and thermal images to derive accurate values of T_c and T_s in sparse canopies remains a challenge.

The forced estimation of H_C using the TSEB-2T_H approach resulted in an overestimation of LE_C and, therefore, hourly transpiration. Given the comparable RMSE obtained from the estimation of T_c using the PT_H and PT_R methods (1.44 and 1.22 °C for PT_H and PT_R), the overestimation deriving from the model could have been attributable to the higher values of R_n,c estimated by the C&N–H model. In contrast to our results, Kool et al. (2021) observed a model tendency to underestimate transpiration due to an underestimation of available energy. To counteract this underestimation, the authors proposed adaptations that included accounting for the higher leaf radiation absorption that would be expected in dense, clumped canopies. In this context, in crops with a complex canopy structure, it is key to apply models capable of accurately estimating canopy transmittance. This is particularly important in productive systems such as the hedgerow system, with high LAI values but low levels of fIPAR. Therefore, one source of the higher error observed in the hedgerow system could be inconsistent clumping index estimates in intensive orchards. Thermal remote sensing models can be particularly sensitive to inhomogeneous distributions of vegetation because clumping affects the relationship between temperature and cover information and the overall energy balance (Kustas and Norman 1997, 1999a, b; Anderson et al. 2005). Many studies have used radiation transfer models and clumping indices in woody row crops (Kustas and Norman 2000; Kustas and Anderson 2009; Parry et al. 2019). But, to our knowledge, these modelling approaches have never been validated in orchards with different production systems and, hence, canopy structures, as is the case of this study. It is possible that the error in transmittance models will be more exacerbated in superintensive production systems with narrower planting distances.

Given the importance in solar irradiance of the PAR spectrum which, according to Campbell and Norman (1998), corresponds to about 45% of the light spectrum, we decided to additionally estimate transpiration with the TSEB-2T_R modeling approach but replacing the C&N-R PAR transmittance model with in situ PAR transmittance measurements (TSEB-2T_fIPAR). Figure 11 shows the regression between sap flow measured transpiration and that obtained with the TSEB-2T_R and TSEB-2T_fIPAR models in trees for which in situ fIPAR measurements were conducted. The results show how the use of in situ fIPAR measurements enhanced transpiration estimates, with the RMSE decreasing from 0.15 mm h⁻¹ (34%) for TSEB-2T_R to 0.10 mm h⁻¹ (29%) for TSEB-2T_fIPAR. Therefore, it occurs to us that the way to improve estimates of transpiration in woody crops in rows can either be to improve the shortwave transmittance model and adapt the clumping index for canopies with different architectures, or directly assimilate the transmittance of the PAR spectrum in the canopy net radiation model.

Fig. 11

Regressions between measured and estimated hourly transpiration with the modelling approaches TSEB-2T_R (a) and TSEB-2T_fIPAR (b). TSEB-2T_fIPAR corresponds to the TSEB-2T_R modelling approach using in situ PAR transmittance measurements. The color indicates each irrigation treatment (TRT), with blue, orange and red representing the fully irrigated (FI), mild stress (MS) and severe stress (SS) treatments, respectively. The point shape shows each production system (PS), with OV, CA and HGR corresponding to open vase, central axis and hedgerow, respectively (color figure online)

Conclusion

As a novel contribution, this study evaluated differences in transpiration estimates in three different production systems and irrigation treatments and undertook validation with sap flow data. In addition, this study compared the performance of transpiration estimation with the TSEB Priestley-Taylor (PT) and contextual (2T) modeling approaches and the effect of applying two clumping indexes to estimate canopy transmittance. One approach considered a basic clumping index for heterogeneous randomly placed clumped canopies (C&N–H) and the other a rectangular hedgerow clumping index (C&N–R).

Significant differences were observed in transpiration between production systems and irrigation treatments. We conclude that the TSEB-2T and C&N-R transmittance models used synergically (TSEB-2T_R) improved transpiration estimations in all production systems and irrigation treatments, with an overall R² of 0.77 and an RMSE of 0.13 mm h⁻¹, representing a relative error of 36%. The better performance of TSEB-2T_R confirms that the availability of high-resolution UAV imagery allows more detailed characterization of the different input parameters needed in the TSEB scheme (f_c, w_c, h_c/w_c). It also allows an adequate separation of canopy and soil temperatures. In addition, the inclusion of estimates of canopy architecture parameters and vegetation indices in a random forest machine learning algorithm allowed a good estimation of LAI. The modeled LAI showed an error of 0.30 m m⁻² in comparison to measured LAI₂₂₀₀.

The better results obtained with the C&N-R shortwave transmittance model highlights the importance of adapting these models to heterogeneous architectures to better estimate canopy transmittance (or R_n,c) and, in consequence, transpiration. Use of a clumping index adjusts the relationship between the overall energy balance, the cover and surface temperature. In our results, one source of error may be an inconsistent estimation of the clumping index in the different almond production systems. For its part, the use of in situ PAR transmittance measurements in the TSEB-2T_R model improved the RMSE from 0.15 to 0.10 mm h⁻¹. Therefore, we suggest that future studies should focus on the improvement and adaptation of shortwave transmittance to heterogeneous canopies and/or on the use of in situ PAR transmittance measurements in canopy net radiation models.