A First Exploration of the Ts/VI “Analytical Triangle” Technique with UAV Imagery for Deriving Key Surface Energy Balance Parameters at Very High Spatial Resolution

Abstract

Knowledge on the spatiotemporal patterns of surface energy balance parameters is crucial for understanding climate system processes. To this end, the assimilation of Earth Observation data with land biosphere models has shown promising results, but they are still hampered by several limitations related to the spatiotemporal resolution of EO sensors and cloud contamination. With the recent developments on Unmanned Aerial Vehicles (UAVs), there is a great opportunity to overcome these challenges and gain knowledge of surface energy balance parameters at unprecedented resolutions. The present study examines, for the first time, the ability of an inversion-modeling scheme, the so-called “analytical triangle” method, to retrieve estimates of surface energy fluxes and soil surface moisture (SSM) at high spatial resolution using UAV data. A further aim of our study was to examine the representativeness of the SSM estimates for the SM measurements taken at different depths. The selected experimental site is an agricultural site of citrus trees located near the city of Palermo on 30 July 2019. The results of comparisons showed that the sensible and latent heat fluxes from UAV were consistent with those measured from the ground, with absolute differences in comparison to ground measurements being 5.00 Wm⁻² for the latent heat (LE) flux and 65.02 Wm⁻² for H flux, whereas for the daytime fluxes H/Rn and LE/Rn were 0.161 and 0.012, respectively. When comparing analytical triangle SSM estimates with SM measurements made at different depths, it was found that there was a gradual increase in underestimation with increasing measurement depth. All in all, this study’s results provide a credible demonstration of the significant potential of the technique investigated herein as a cost-effective and rapid solution for estimating key parameters characterizing land surface processes. As those parameters are required by a wide range of disciplines and applications, utilization of the investigated technique in research and practical applications is expected to be seen in the future.

1. Introduction

Earth’s surface is characterized by the constant energy and mass interactions occurring in the atmosphere-biosphere system [1,2,3,4]. Climate change is a factor that raises several issues that are either local or supra-local in nature. These issues are related with soil degradation, resource (water) scarcity, food scarcity observed in many parts of the world, etc. It is therefore imperative to have a deeper knowledge of all of these natural processes occurring on the planet [5,6]. They include the physical processes of highly spatiotemporal variables such as surface-soil moisture (SSM) and latent (LE) and sensible heat (H) fluxes. These variables either directly or indirectly affect the water cycle and, therefore, the flow of ecosystems with regard to their physical processes [7,8,9,10]. One such influence, for example, is the way in which the soil-atmosphere system regulates the partitioning of the energy present at the soil upper layer into LE and H through evaporation and transpiration processes, thereby linking the water-energy balance to the soil’s moisture–temperature dipole [11,12]. The latter two variables are very important: H flux alters the atmospheric turbulence at the Earth’s surface, whereas LE flux affects the global water and carbon cycle [13,14]. Therefore, their accurate estimation is of prime interest for a number of environmental and commercial applications, from sustainable water resource management to evaluating parameterization schemes for weather and climatic models [15,16,17].

The advent of Earth Observation (EO) has provided economically feasible means to derive temporally consistent coverage of those parameters at different geographical scales [18,19]. To this end, Unmanned Aerial Vehicles (UAVs) have attracted the attention of researchers for monitoring energy surface parameters [20,21,22]. Compared to satellites, UAVs are much more adaptive, agile and versatile, and they could offer rapid and reproducible deployment for collecting high-resolution data near real-time with low operational costs [23,24,25]. Due to their lower operating costs, near real-time image acquisition and increased spatiotemporal resolution, UAVs are a desirable system for the cartographic representation and continual surveillance of energy and water surface parameters ([26,27,28]). UAVs offer spatial resolutions ranging at sub-meter, meter, and centimeter levels, dependent on the altitude of flight. This variability is particularly advantageous for assessing cultivated plants like vines and trees with clustered canopy structures [29].

Apart from the increased availability in suitable datasets and platforms that are being used for monitoring energy surface parameters, several methodologies have been developed, varying from statistical-semi-empirical to fully analytical ones that simulate physical processes and are fully physics-based [30]. These modeling schemes are characterized by their varying mechanisms and degrees of complexity, data requirements, basic assumptions, and accuracy. Evidently, there is a specific group of EO-based techniques that aim to deduce the surface fluxes of LE, H, and/or SSM at a variety of spatial and temporal scales based on the synergy of satellite data from optical (visible and infrared—VNIR) and thermal infrared (TIR) radiometers These methods are known as Ts/VI and they essentially simulate the physical relationship that land temperature (Ts) has with a spectral index related to vegetation (VI) and associated with fractional vegetation cover (Fr) [31,32,33].

A number of research studies published in the scientific literature have demonstrated the existence of a triangular (or trapezoidal) shape of the Ts/VI space that emerges from a scatterplot [34,35,36,37]. Such a scatter plot shows a triangular (or trapezoidal) shape that has four distinct physical boundaries. The “dry edge” or “warm edge” is defined by the points of highest temperature that contain differing amounts of vegetation and bare soil. It is assumed to represent conditions of limited surface-soil water content and zero evaporative flux from the soil, characterizing surfaces with the largest water stress for a range of VIs. Likewise, the “cold” or “wet” edge depicts those areas where they have the highest soil moisture and therefore potential SSM/ET for a range of VIs [38]. For pixels with the same VI, those with minimum Ts represent the case of the strongest evaporative cooling, whereas those with maximum Ts represent those with the weakest evaporative cooling [39]. Detailed descriptions on the information encapsulated in the Ts/VI triangle domain physical properties as well as the main factors driving the shape of the Ts/VI scatterplot are summarized in [40].

It has already been established that the retrieval of such spatially distributed fluxes and/or SSM using the Ts/VI ‘triangular’ scatterplot is achievable without the use of a boundary-layer model. Yet, more sophisticated approaches involve a land biosphere model, specifically a Soil Vegetation Atmosphere Transfer (SVAT) model, via a technique commonly termed as the “analytical triangle” [34]. Such approaches aim to combine the vertical coverage and fine temporal continuity of those models with horizontal coverage and spectral resolution of EO data. Various validation studies of the “analytical triangle” technique have demonstrated its ability to provide estimates of both surface heat fluxes and SSM with accuracies of 40–70 Wm⁻² and within 5% vol vol⁻¹ for SSM over homogenous areas [36,39,41]. Yet, validation exercises on this method are scarce and so far most of them have been implemented using either satellite data from coarse to medium spatial resolution sensors or from medium to relatively high spatial resolution airborne sensing systems.

Also, notably, most of the technique validation exercises have so far focused on comparing the predicted SSM against the surface-layer-soil moisture from reference data coming mainly from in situ measurements. Hence very little is known at present on the agreement between the “triangle”-predicted SSM and reference measurements at variant soil depths. As has been noted in several studies, SSM retrievals based on the triangular feature space are sensitive to a very shallow soil layer (e.g., [42]). All of the above points raise a high degree of interest for further studies on the “analytical triangle” method, which is also very topical given that variants of this method are currently being explored for operational implementation by different agencies.

Hence, in view of the considerations above, the present study’s objectives are two-fold: (i) to study, for the first time, the use of the “analytical triangle” method for the retrieval of turbulent energy flux parameters using very high-resolution images acquired by a drone in a Mediterranean environment, and (ii) to assess the agreement of the “analytical triangle” SSM estimates with collocated in situ measurements taken at different depths.

2. Experimental Set Up

2.1. Study Site and Ground Instrumentation

Our selected experimental site is in the vicinity of Palermo city, Italy (38°4′53.4″ N, 13°25′8.2″ E). The site covers an area of approximately one hectare with very small variation in altitude, ranging from 25–50 m above sea level, and with slopes from 1–3%. The site is representative of an eastern Mediterranean region, with the primary land use being citrus trees with understory grass seasonally covering the bare soil.

This site is owned and managed by the University of Palermo and is equipped with very detailed ground instrumentation that continuously monitors several environmental and biophysical variables. At the site has been installed an eddy covariance system for the measurement of radiation and turbulent fluxes. The eddy covariance InfraRed Gas Analyzer and a Sonic Anemometer have been installed slightly above the vegetation canopy, at 3.5 m above ground level (a.g.l.), i.e., approximately at 55 and 95 cm above the vegetation canopy. Furthermore, a meteorological station is also installed that records multiple micrometeorological parameters such as relative wind speed and direction, rainfall, solar radiation, air humidity, and temperature. A four-component net radiometer is used to measure the individual components of net radiation (R_n). In addition, a soil-moisture measurement network is deployed at the site (namely drill and drop frequency domain reflectometry sensors), in a total of eight sites to monitor soil moisture every 10-cm from the top 5 cm of to 60 cm of the soil profile. As can be observed in Figure 1, four of the sites, specifically Sites 1 to 4, are installed over well-irrigated areas, and Sites 5–8 are installed over a deficit-irrigation scheme. To facilitate data collection, a communication station enabling remote access to the sensors network datasets has been installed at the experimental site.

2.2. UAV Data: Acquisition

UAV imagery was acquired within the bounds of the experimental site, on 30 July 2019, at around 1200 noon (GMT). The dataset consisted of nine multispectral bands with a 10 nm bandwidth, which were acquired with a 36.5° field of view using an NT8 contras octocopter carrying a RikolaDT-17 camera (Rikola Ltd., Oulu, Finland). Thermal imagery was acquired using a DJI Mavic 2 Enterprise Dual quadcopter carrying on-board a FLIR Lepton^® (FLIR Systems, Inc., Wilsonville, OR, USA) acquiring the data in the longwave infrared spectral range (from 8 to 14 µm). The average GSD was 3.46 cm. Nine black-and-white control targets were placed throughout the study area and used in a Global Navigation Satellite System (GNSS) survey. The coordinates of the targets were measured using a Topcon Hiper V receiver and Network Real-Time Kinematic (NRTK) positioning, available by the UNIPA GNSS CORS network. Furthermore, field-based spectroradiometric measurements were collected at sample locations composed of both bare soil and vegetation and a handheld thermograph was also employed to measure surface temperature independently of the UAV data acquisitions. In situ surface temperature measurements (Ts) were carried out at noon using a handheld FLIR SC660. Further details concerning the in situ and UAV data acquisition at the experimental site are available, for example in [39].

3. Methods

The methodological approach implemented in this study can be summarized in the following three stages: (i) UAV data pre-processing, including the estimation of fractional vegetation cover (Fr) and normalization of the UAV-derived temperature datasets, (ii) implementation of the triangle technique, and (iii) accuracy assessment of the predictions made by the investigated technique. A graphical illustration of the methodological approach employed to satisfy this study’s objectives is illustrated in Figure 2. In the following sections, more details are provided concerning the implementation of the presented herein methodological approach methodology.

3.1. UAV Data Pre-Processing

Standard pre-processing included the data calibration to radiance or reflectance, georeferencing, bands co-registration, and image resampling to match the spatial resolution between optical and TIR images. These include the UAV image data orthorectification and their resampling at 4 cm spatial resolution. The resampling scale has been based on the average GSD (see Section 2.2). Regarding the orthorectification, a standard photogrammetric/SfM approach was utilized, using Pix4D mapper (by Pix4D Inc.). As a next step, the empirical line technique was implemented in the VIS-NIR imagery to convert digital numbers to ground reflectance [43]. Similarly, a linear regression with at-ground thermographs and an emissivity map of the soil vegetation system were used to calibrate TIR images into surface radiometric temperatures [44]. Emissivity was computed by adopting the approach proposed by Valor and Caselles [45], while assuming emissivity values of bare soil and densely vegetated ground to be equal to 0.97 and 0.99, respectively [46]. More details on these pre-processing steps can be found in [39].

Following this step, Fr was estimated from the Normalized Difference Vegetation Index (NDVI). This was computed by first scaling NDVI to an N* value calculated using the following formula:

$${\mathrm{N}}^{*}=\frac{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}-{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}_0}{{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}_{\mathrm{s}}-{\mathrm{N}\mathrm{D}\mathrm{V}\mathrm{I}}_0}$$

(1)

where NDVI_s and NDVI₀ are the maximum and minimum NDVI values at minimum (0%) and maximum (100%) vegetation covers, respectively.

Then, the transformation of scaled NDVI (N*) to Fr is calculated as follows:

$$\mathrm{F}\mathrm{r}={\mathrm{N}}^{*2}$$

(2)

where N* is the linearly scaled NDVI and Fr is the vegetation fraction. Transformation of N* to Fr allows the measured surface radiant temperatures from the UAV data products and the SVAT-simulated to be plotted on the same scale.

After the vegetation normalization, the next step included the temperature normalization using the following formula (Equation (3)):

$${\mathrm{T}}_{\mathrm{s}\mathrm{c}\mathrm{a}\mathrm{l}\mathrm{e}\mathrm{d}}=\frac{{\mathrm{T}}_{\mathrm{o}}-{\mathrm{T}}_{\mathrm{m}\mathrm{i}\mathrm{n}}}{{\mathrm{T}}_{\mathrm{m}\mathrm{a}\mathrm{x}}-{\mathrm{T}}_{\mathrm{m}\mathrm{i}\mathrm{n}}}$$

(3)

where T_max and T_min are the expected maximum and the minimum T for wet vegetated soil and for bare-soil pixels, respectively, as interpolated from a scatterplot of UAV-based derived T versus Fr. For any pixel in the scene, T_o corresponds to the T value. Τhe degree to which the surface temperature resembles that of either a wet, vegetated surface or a dry, bare-soil surface is expressed by T_scaled, which is a generalized, unitless, measure of temperature. T_scaled varies with time of day and atmospheric conditions [34]. Therefore, T_scaled is a more fundamental variable than Ts is.

The final UAV image products derived upon completion of the pre-processing steps described above, namely the T_scaled and Fr maps, are illustrated in Figure 3.

3.2. Analytical Triangle Implementation

SimSphere is the SVAT model that was used for the analytical triangle implementation. Briefly, SimSphere is a 1-dimensional SVAT model with a plant component [47,48]. It can simulate several land–atmosphere exchanges of energy, water and other gases that are taking place in a column that extends from the root zone below the soil surface up to a level well above the surface canopy, which is the top of the surface mixing layer. In the model, soil fluxes and vegetation meld at the top of the vegetation canopy. Their relative weights depend on Fr, which is defined by the user. SimSphere simulates land–atmosphere processes over a 24-h cycle at a chosen time step, starting from a set of initial conditions given in the early morning (at 06:00 a.m. local time) with continuously evolving interactions between soil, plant, and atmosphere layers. The model parameterization requires over 50 input parameters in total, which are categorized into 7 defined groups: location and time, surface, vegetation, meteorological, hydrological, atmospheric and soil. SimSphere assesses the diurnal evolution of more than 30 prognostic variables associated with the radiative, hydrological, and atmospheric physical domains. Several other outputs are also predicted such as the wind velocity, air temperature, humidity, and radiometric surface temperature at various levels in and above the canopy, plus several plant parameters such as leaf water potential and stomatal resistance. Detailed insights into the recent architectural advancements performed in the model is provided in [49]. The latest model version can be accessed at no cost through the Department of Meteorology at Pennsylvania State University, USA (https://simsphere.ems.psu.edu/model.html (accessed on 19 June 2024)).

A comprehensive account of the analytical triangle method is provided in [34]. Briefly, the method allows the retrieval of the LE and H fluxes and soil-water availability (Mo), which represents surface wetness in the bare-soil surface (top few millimeters of it). To derive the inversion equations required to estimate the spatially explicit maps of land surface LE/H fluxes and Mo by the analytical triangle method, the method requires the use of a SVAT model, which is coupled to the UAV data. In this study, the latest version of the SimSphere SVAT model was used for this purpose. SimSphere was parameterized using the time and geographical location and the general soil and vegetation characteristics of our experimental site, together with the appropriate atmospheric profile data acquired at no cost from the University of Wyoming database, (http://weather.uwyo.edu/upperair/sounding.html accessed on 19 June 2024). Parameters were adjusted based principally on the available ground measurements and ancillary information that was obtained directly from the site manager (pers. comm. Professor Provenzano). Next, SimSphere was iterated repetitively until the extreme values of Fr and Ts in the Ts/Fr scatterplot between the simulated (from SimSphere) and observed (UAV-derived data) were matched (i.e., initial model simulations endeavored to align observed Ts with two endpoints (NDVIo, NDVIs) where they intersect the with “dry” edge). This extrapolation of NDVIo and NDVIs guarantees that the implied temperatures along the “dry” edge for bare soil and full Fr are consistent with simulations for an Mo of zero.

Once SimSphere tuning was completed for the test-site image’s date, the model was run several times while kee** the time (corresponding to UAV overpass) the same but varying Fr and Mo over all possible values (0–100% and 0–1, respectively) in increments of 10 and 0.1 respectively. The result was a matrix of model outputs, comprising Mo, Fr, T_scaled, LE, and H, for the time of the UAV image acquisition (i.e., at noon UTC time) that were calculated for each combination of Fr and Mo. Next, the output matrix computed from the previous step was used to derive a series of non-linear (quadratic) equations, empirically relating Fr and Tscaled to each of the other variables of interest: H, LE, LE/Rn and H/Rn. Through this process, the set of physically based relationships between the various surface-atmosphere parameters, as described by the detailed biophysical descriptions included in SimSphere and inherent in the matrix outputs, were used to derive a series of simple, empirical relations relating each of these parameters to only the Fr and T_scaled values recorded at that particular location. As these variables of Fr and T_scaled are derivable from the UAV data, these empirical equations are then used to derive the required spatially explicit maps of the land surfaces’ LE and H fluxes as well as Mo from the UAV-derived Fr and T_scaled values.

The quadratic polynomial equations derived from the SVAT model matrix model outputs have a general form as follows (shown here for the version relating M₀ to Fr and T_scaled):

$${\mathrm{M}}_0={\displaystyle \sum _{\mathrm{p}=0}^3}{\displaystyle \sum _{\mathrm{q}=0}^3{\mathrm{a}}_{\mathrm{pq}}{({\mathrm{T}}_{\mathrm{scaled}})}^{\mathrm{p}}{({\mathrm{F}}_{\mathrm{r}})}^{\mathrm{q}}}$$

(4)

where the coefficients p and q are derived from non-linear regression between the matrix values of Fr, T_scaled, and Mo as they vary from 0 to 3.

In this matrix, polynomial equations were calculated to empirically relate the following variables to the values of Fr, Tscaled, LE, H, LE/Rn, and H/Rn. The latter two parameters are reportedly close approximations of the daily average expressions as they are less dependent of the time of the day (because they remain relatively constant throughout the daytime period under cloud-free conditions) and are related to Mo and Fr [41,50]. The last step of the method’s implementation involved converting Mo values to SSM values by multiplying Mo by the soil’s field capacity (this information was obtained by the site managers).

3.3. Accuracy Assessment

The availability of the ground measurements that were almost concurrently collected at the experimental site during UAV data acquisition allowed the quantitative evaluation of the “analytical triangle” outputs. At the experimental site, ground measurements of the radiation and turbulent fluxes were available for one location only, whereas soil-moisture measurements were available from a total of eight sites across the field on which different irrigation treatments were applied such that it was possible to compute a series of statistical metrics. Therefore, for SSM, the agreement between predictions and observations were assessed based on a series of appropriate statistical metrics (

Table 1. The statistical metrics employed herein are used to validate the analytical triangle retrievals and the ground reference. Subscripts i = 1 … N denote the individual observations, P denotes the predicted values, and O denotes the “observed” values, which in our case are those obtained from our experimental site. The horizontal bar in the scatter ratio equation denotes the mean value.

Name	Description	Mathematical Definition
Bias	Bias (accuracy)	$\mathrm{b}\mathrm{i}\mathrm{a}\mathrm{s}=\frac{1}{\mathrm{N}}{\displaystyle \sum _{\mathrm{i}=1}^{\mathrm{N}}}({\mathrm{P}}_{\mathrm{i}}-{\mathrm{O}}_{\mathrm{i}})$
Scatter	Scatter (precision)	$\mathrm{s}\mathrm{c}\mathrm{a}\mathrm{t}\mathrm{t}\mathrm{e}\mathrm{r}=\frac{1}{(\mathrm{N}-1)}{\displaystyle \sum _{\mathrm{i}=1}^{\mathrm{N}}}\sqrt{({\mathrm{P}}_{\mathrm{i}}-{\mathrm{O}}_{\mathrm{i}}\overline{({\mathrm{P}}_{\mathrm{i}}-{\mathrm{O}}_{\mathrm{i}})}{)}^2}$
RMSD	Root Mean Square Difference	$\mathrm{R}\mathrm{M}\mathrm{S}\mathrm{D}=\sqrt{\mathrm{b}\mathrm{i}\mathrm{a}{\mathrm{s}}^2+\mathrm{s}\mathrm{c}\mathrm{a}\mathrm{t}\mathrm{t}\mathrm{e}{\mathrm{r}}^2}$

). Those included the root mean square difference (RMSD), the scatter. and the bias. These statistical metrics have been widely and prominently applied in analogous validation studies from satellite data, including previous verification exercises of the “analytical triangle” technique (e.g., [36,38,39,41,51]).

4. Results

The map products derived from the “analytical triangle” implementation and their corresponding histograms are illustrated in Figure 4. These products were included the instantaneous Mo, LE, and H fluxes as well as in the daytime average fluxes of LE/R_n and H/R_n. A first visual inspection of the derived maps (Figure 4) shows a reasonable range of values for the different outputs within the experimental site. The spatial patterns of the derived maps are generally arranged by the existing land cover classes and the Fr and Ts UAV-based maps of the experimental site. For example, it is observed that LE, SSM, and LE/R_n values are lower in the non-vegetated surfaces than in the vegetated surfaces, whereas the observed distributions of H and H/R_n values are the opposite. Furthermore, the predicted parameters show a reasonable amount of spatial variability between them (e.g., when comparing M_o against LE and H flux measures), over both the vegetated and the bare-soil areas.

Table 2. Summary of the point-by-point comparisons conducted.

Flux Comparisons	Predicted (P)	Observed (O)	Difference (P-O)
LE flux (Wm−2)	193.62	198.62	−5.00
H flux (Wm−2)	271.67	336.69	−65.02
LE/Rn	0.325	0.313	0.012
H/Rn	0.456	0.295	0.161

summarizes the results from the point comparisons that were performed for the different variables predicted by the triangle technique. As can be observed from Table 2, the “analytical triangle” provided satisfactory predictions of other fluxes, as the absolute differences compared to ground measurements were 5.00 Wm⁻² for LE flux, 65.02 Wm⁻² for H flux, 0.161 for H/R_n, and 0.012 for LE/R_n. In the overall evaluation of the results, spatial differences between the “analytical triangle” predictions and the ground reference data can cause large mismatches when compared and, thus, influence the produced results. Uncertainties or errors thus arise from the lack of representativeness in data due to scale mismatch. Finally, uncertainty in the in situ observations themselves should be considered when evaluating the results.

On the other hand, to assess the accuracy of the SSM predictions by the “triangle”, in situ soil-moisture measurements were collected at six different depths across eight different locations on the experimental site. The statistical analysis results obtained from the direct comparison are summarized in Figure 5 below.

The comprehensive comparison of the analytical triangle SSM estimates with the in situ SSM measurements taken from 10 to 60 cm has shown a slight overestimation in the upper layer (10 cm) values with a positive bias of 0.004 cm³ cm⁻³. However, a decrease in the bias is observed, indicating a gradual rise in the underestimation, as the measurement depth increases, reaching a bias value of −0.270 cm³ cm−³ in the measurement depth of 60 cm, which is possibly influenced by factors such as soil composition and the measurement. This trend is similarly observed in RMSD, where the lowest value of RMSD (0.092 cm³ cm⁻³) is observed in the upper soil depth that gradually increases as the soil depth increases, finally reaching a value of 0.302 cm³ cm⁻³ RMSD at the lowest measurement depth (60 cm). In terms of scatter values, the lowest scatter was observed at the 40 cm measurement depth with a value of 0.087 cm³ cm⁻³, whereas the highest scatter value was observed, similarly with bias and RMSD, at the measurement depth of 60 cm.

5. Discussion

The present study investigated, to our knowledge for the first time, the performance of the so-called “analytical” triangle when used with very high-resolution UAV imagery, which was acquired at an arid/semi-arid environment located in Palermo, Italy. In this framework, it also investigated the agreement between the SSM predictions made by the triangle with ground measurements obtained at different soil depths.

Overall, the results obtained from the present study, even though they are based on single-image analysis, confirm the usefulness of the examined technique for SSM/M_o and energy flux retrievals at very fine spatial resolution when implemented with UAV data. The spatial patterns of the derived maps are in accordance with the different landscape characteristics of the study area such as land cover type and topography. Statistical that was analysis also performed using the available ground measurements showed that the “analytical triangle” was able to satisfactorily predict all of the derived parameters. Notably, statistical analysis results reported herein are in close agreement with those reported by other studies retrieving SSM and LE/H flux values using TIR-based techniques (e.g., [36,38,39]. In addition, particularly for the SSM comparisons at the different depths, after the M_o values were converted to SSM for the eight stations, the results showed that accuracy was gradually reduced as the SSM measurement depth increased. Those results are also in accordance with findings from similar studies and further confirm that SSM estimates from the analytical triangle implementation are representative of the upper soil layer [42,48]. Furthermore, our study results are also in accordance with other studies utilizing UAV datasets to monitor energy surface parameters. Ref. [23] utilized the Two-Source Energy Balance (TSEB) land surface scheme with UAV imagery and reported similar accuracies in H flux (−58 W m⁻²) but with noticeable differences in LE flux (63 W m⁻²).

In the statistical results reported herein, spatial discrepancies between the “analytical triangle” estimates and the ground reference values can cause large mismatches when compared. Uncertainty in the in situ observations themselves should be considered when evaluating the results. These may stream from the lack of representativeness in the reference data used to validate the predictions due to scale-mismatch errors, geo-location errors, and noise from surface heterogeneity. In addition, the lack of complete agreement might be in part related to uncertainty in the retrieval of Fr and Ts values from the UAV data [39]. In addition, uncertainties are also included in the conversion of M_o to SSM values, supposing a constant value of field capacity (equal to 0.24) for the entire field. Field capacity is a highly dynamic variable both in spatial and temporal aspects. Obtaining knowledge of the spatial distribution of field capacity would improve the conversion of Mo to SSM, thus reducing the uncertainty in the comparison of the two datasets.

One of the key advantages of the “analytical triangle” includes the use of the SVAT model, which allows combining of the horizontal coverage and spectral resolution of EO data with the vertical coverage and fine temporal continuity of the model. Such approaches combine the horizontal coverage and spectral resolution of EO data with the vertical coverage and fine temporal continuity of those models [49]. Also, the “analytical triangle”, in contrast to other Ts/VI techniques, provides a non-linear interpretation of the Ts/VI space and, thus, a solution for the computation of the energy flux values and/or SSM, which can be a more realistic assumption of flux and SSM dynamics. What is more, it offers a relatively easy transformation of the instantaneous-derived energy fluxes to expressions of daytime averages [40].

Yet, further work is required to help establish the full potential of the technique. This involves, among other things, its extensive implementation in different ecosystems or the examination of the atmospheric sounding effect on the inversion results. Further research will be conducted to explore the effects of using different methods for Fr estimation from the UAV data on the variables predicted by the methods examined herein. This would provide information on the sensitivity of the method to the Fr input parameter.

6. Conclusions

In summary the present study assesses the performance of the so-called “analytical triangle” technique when applied with very high spatial resolution UAV imagery, using as a test site a typical Mediterranean orchard-field environment located in Sicily, Italy. A further objective was to specifically evaluate the agreement of the SSM predicted by the “analytical triangle” against concurrent SSM ground measurements available at different depths at our experimental site. The innovation introduced in the present study is highly significant as, to our knowledge, our study represents the first effort in applying this specific inversion modeling with UAV data while also incorporating verification of the soil-depth-analysis components.

All in all, the results demonstrated good performance of the investigated technique in reproducing turbulent-daytime surface energy parameters for the study area, and in creating SSM maps at high spatial resolution and with accuracy comparable to, or better than, that reported in analogous verification studies that used low- to medium-resolution satellite imagery. Furthermore, the results confirmed for the first time with such high-resolution data that the SSM estimates are representative of upper-layer-soil moisture. In conclusion, this study credibly demonstrated the significant potential of the “analytical triangle” method in deriving high-resolution estimates of LE/H flux and SSM utilizing high-resolution UAV data. The combination of UAV imagery with this technique provides a cost-effective and rapid solution for estimating key parameters characterizing land surface processes, which are required in a wide range of disciplines and applications. Our study also provides further evidence supportive of placing the “analytical triangle” technique in a privileged position as a candidate for further examination for prospective operationalization with either satellite data or airborne-EO data.