1.2.3. Fractional Cover (*fc*) and Canopy Width (*wc*):

1.2.3. Fractional Cover (*fc*) and Canopy Width (*wc*): Fractional cover (*fc*) is the proportional area of vine for each spatial domain under analysis, where values vary from 0 through 1. *fc* is used to estimate *wc* and the clumping index, which is a factor to adjust the remotely sensed *LAI* value, which is assumed to be uniformly distributed (homogeneous) over the landscape instead of being clumped [54]. These are used to estimate the actual canopy gap fraction, which is greater than the homogenous case. It is required as input for the radiation transmission and wind extinction algorithms through the canopy layer. The magnitude of *wc* is a length scale representing the area occupied by vine leaves along the vine row, which varies spatially and temporally based on phenology and management (i.e., vine manipulation via the trellis system and pruning) (Figure 2). Fractional cover (*fc*) is the proportional area of vine for each spatial domain under analysis, where values vary from 0 through 1. *f<sup>c</sup>* is used to estimate *w<sup>c</sup>* and the clumping index, which is a factor to adjust the remotely sensed *LAI* value, which is assumed to be uniformly distributed (homogeneous) over the landscape instead of being clumped [54]. These are used to estimate the actual canopy gap fraction, which is greater than the homogenous case. It is required as input for the radiation transmission and wind extinction algorithms through the canopy layer. The magnitude of *w<sup>c</sup>* is a length scale representing the area occupied by vine leaves along the vine row, which varies spatially and temporally based on phenology and management (i.e., vine manipulation via the trellis system and pruning) (Figure 2).

### 1.2.4. *wc*/*h<sup>c</sup>* Ratio

1.2.4. *wc/hc* Ratio In *TSEB* and *TSEB2T* models, the *wc*/*hc* ratio is required as input to the radiation transmission and wind extinction algorithms through the canopy layer developed for vineyards [2,55]. The *wc*/*hc* In *TSEB* and *TSEB2T* models, the *wc*/*h<sup>c</sup>* ratio is required as input to the radiation transmission and wind extinction algorithms through the canopy layer developed for vineyards [2,55]. The *wc*/*h<sup>c</sup>* ratio value is obtained by simply calculating canopy width over canopy height (Figure 2).

ratio value is obtained by simply calculating canopy width over canopy height (Figure 2).

growing seasons (May–August) using different spatial domain scales.

#### **2. Materials and Methods**

**2. Materials and Methods**  The methodology to assess the impact of changes in the contextual spatial domain for the *TSEB2T* model is graphically presented in Figure 3. The analysis was performed for wine grape The methodology to assess the impact of changes in the contextual spatial domain for the *TSEB2T* model is graphically presented in Figure 3. The analysis was performed for wine grape growing seasons (May–August) using different spatial domain scales.

**Figure 3.** Study methodology for assessing the impact of the *TSEB2T* contextual spatial domain. **Figure 3.** Study methodology for assessing the impact of the *TSEB2T* contextual spatial domain.

#### *2.1. Study Area and Data Sources*

*2.1. Study Area and Data Sources*  The study site is located near Lodi, California (38.29୭N, 121.12୭W) with an area of approximately 150 ha. The two vineyard blocks (north and south) are part of the Sierra Loma vineyard ranch (Figure 4). The north block was planted in 2009, while the south block was implemented in 2011, leading to different levels of vine maturity, and hence, biomass and grape production. Both vineyards are managed cooperatively by Pacific Agri-Lands Management. The plantation structure in both fields is the same, with vine rows having east–west orientation with a row width of 3.35 m (11 feet). A cover crop grows in the interrows, occupying ~2 m, with bare soil strips along the vine rows spanning ~0.7 m. The purpose of the cover crop is to deplete plant available water in the interrows from the fall and winter precipitation in order to control vine growth in the spring by irrigation. Typically, the vine height varies between 2 m and 2.5 m above ground level (*agl*) and vine biomass is concentrated mainly in the upper half of the vine canopy height. The actual vine canopy width varies spatially and temporally due to vine management practices. This study site is a part of the Grape Remote Sensing Atmospheric Profile and Evapotranspiration eXperiment (*GRAPEX*) project run by the USDA Agricultural Research Service in collaboration with E&J Gallo Winery, Utah State University, University of California in Davis, and others [56]. The study site is located near Lodi, California (38.29◦N, 121.12◦W) with an area of approximately 150 ha. The two vineyard blocks (north and south) are part of the Sierra Loma vineyard ranch (Figure 4). The north block was planted in 2009, while the south block was implemented in 2011, leading to different levels of vine maturity, and hence, biomass and grape production. Both vineyards are managed cooperatively by Pacific Agri-Lands Management. The plantation structure in both fields is the same, with vine rows having east–west orientation with a row width of 3.35 m (11 feet). A cover crop grows in the interrows, occupying ~2 m, with bare soil strips along the vine rows spanning ~0.7 m. The purpose of the cover crop is to deplete plant available water in the interrows from the fall and winter precipitation in order to control vine growth in the spring by irrigation. Typically, the vine height varies between 2 m and 2.5 m above ground level (*agl*) and vine biomass is concentrated mainly in the upper half of the vine canopy height. The actual vine canopy width varies spatially and temporally due to vine management practices. This study site is a part of the Grape Remote Sensing Atmospheric Profile and Evapotranspiration eXperiment (*GRAPEX*) project run by the USDA Agricultural Research Service in collaboration with E&J Gallo Winery, Utah State University, University of California in Davis, and others [56].

calibrated [58].

are provided by Alfieri et al. and Agam et al. [59,60].

**Figure 4.** Layout of study area in Lodi, California, locations of *EC* towers and example of 90% of *EC* footprint at afternoon for 02 June 2015. **Figure 4.** Layout of study area in Lodi, California, locations of *EC* towers and example of 90% of *EC* footprint at afternoon for 02 June 2015.

Flights campaigns were conducted by the *AggieAir sUAS* program at Utah State University (https://uwrl.usu.edu/aggieair/). Optical and thermal high-resolution imagery of the study site were collected from different flights in 2014, 2015, and 2016. Vegetative and soil conditions changed between the field campaigns. The 2016 flight imagery represents the early part of the growing season, around the time phenologically of fruit set, while other flights in 2014 and 2015 represent full vine canopy development and grape vine phenology in the pre- and post-veraison stages. Table 1 lists information concerning the different flights. The pixel resolution of the *sUAS* imagery collected is 10 cm and 60 cm for the optical and thermal bands, respectively. The spectral range of the optical data is similar to Landsat and includes visible bands (red, green, and blue) as well as near-infrared. However, the thermal band is different than Landsat, having a bandwidth spanning from 7 to 14 µm [57]. Thermal data, acquired using a lightweight micro-bolometer camera, were radiometrically Flights campaigns were conducted by the *AggieAir sUAS* program at Utah State University (https://uwrl.usu.edu/aggieair/). Optical and thermal high-resolution imagery of the study site were collected from different flights in 2014, 2015, and 2016. Vegetative and soil conditions changed between the field campaigns. The 2016 flight imagery represents the early part of the growing season, around the time phenologically of fruit set, while other flights in 2014 and 2015 represent full vine canopy development and grape vine phenology in the pre- and post-veraison stages. Table 1 lists information concerning the different flights. The pixel resolution of the *sUAS* imagery collected is 10 cm and 60 cm for the optical and thermal bands, respectively. The spectral range of the optical data is similar to Landsat and includes visible bands (red, green, and blue) as well as near-infrared. However, the thermal band is different than Landsat, having a bandwidth spanning from 7 to 14 µm [57]. Thermal data, acquired using a lightweight micro-bolometer camera, were radiometrically calibrated [58].

To evaluate the *ET* performance at different spatial domain scales, two eddy covariance (*EC*) flux systems were deployed for the measurements of turbulent fluxes, including *LE* and *H*, and the available energy terms of *Rn* and *G*. Both towers are located at the eastern edge of the fields, due to predominant winds from the west. Ground measurements, including soil temperature and soil moisture were also collected. A complete listing of all measurements on the towers is given by Kustas et al. [56]. Details of the post processing of the *EC* data as well as the available energy measurements


**Table 1.** Dates and times of AggieAir *GRAPEX* flights used in this study.

To evaluate the *ET* performance at different spatial domain scales, two eddy covariance (*EC*) flux systems were deployed for the measurements of turbulent fluxes, including *LE* and *H*, and the available energy terms of *R<sup>n</sup>* and *G*. Both towers are located at the eastern edge of the fields, due to predominant winds from the west. Ground measurements, including soil temperature and soil moisture were also collected. A complete listing of all measurements on the towers is given by Kustas et al. [56]. Details of the post processing of the *EC* data as well as the available energy measurements are provided by Alfieri et al. and Agam et al. [59,60].

*EC* micrometeorological data also included wind speed, air temperature, vapor pressure, air pressure, and shortwave radiation. Hourly average values of these atmospheric forcing variables, as well as the components of the surface energy balance, were computed. Table 2 illustrates the in-situ micrometeorological parameters and the name of the instruments used for the measurements.


**Table 2.** Description of in-situ micrometeorological measurements in this study.

<sup>1</sup> The use of trade, firm, or corporation names in this article is for the information and convenience of the reader. Such use does not constitute official endorsement or approval by the US Department of Agriculture or the Agricultural Research Service of any product or service to the exclusion of others that may be suitable.

Given the high fluctuation of atmospheric conditions during the daytime, the flux footprint or contributing source area of each *EC* tower was estimated for the hourly period encompassing *sUAS* flight campaigns using the two-dimensional (*2D*) flux footprint model developed recently by Kljun et al. [61]. Because a 100% *EC* footprint fetch could extend over the study area, a 90% footprint area (90% cutoff) was used for analysis. Then, the weighted footprint area was divided by 0.9.

#### *2.2. Data Processing*

In this study, images were acquired remotely by *sUAS*, and the data were terrain corrected using georeferencing based on ground control points (*GCPs*). Furthermore, both thermal and optical data were atmospherically corrected.

#### 2.2.1. Thermal Data

Torres-Rua [57] indicated that the thermal data obtained from the *sUAS* thermal sensors in this study are adversely affected by changes in transmissivity and atmospheric radiance. For this reason, ground measurements of temperature were collected in the same timeframe as the *sUAS* flight and compared with the imagery to calibrate the thermal image data. More details about the calibration of temperature imagery related to this study can be found in Torres-Rua [57].

#### 2.2.2. Optical Data

Radiometric agreement between remotely sensed data from different platforms constitutes one of the major challenges in image processing. Therefore, in this research, the images acquired by *sUAS* were upscaled and harmonized with Landsat using the point spread function (*PSF*). More details related to *sUAS* data harmonization can be found in Hassan-Esfahani et al. [62].

#### *2.3. Energy Balance Closure Adjustment Methods for EC*

While the *EC* technique provides measurements of turbulent fluxes *H* and *LE*, a lack of energy balance closure with the available energy terms *R<sup>n</sup>* and *G* [63] is well documented. This results in *R<sup>n</sup>* − *G* > *LE* + *H* [64,65], and the computed closure ratio (*CR*) evaluates the energy balance discrepancy, *CR* = (*H* + *LE*)/(*R<sup>n</sup>* − *G*). This ratio varies during the daytime, but for the *sUAS* flights [55] it was found to be above 0.8, except for the May 2 afternoon flight where it fell to around 0.7.

To avoid any bias when comparing the energy balance models with *EC* field measurements, the energy closure issue needs to be handled and resolved. Twine et al. [66] suggested a method for energy balance closure that assumes the Bowen ratio (*H*/*LE*) before and after adjustment are the same, while considering both *R<sup>n</sup>* and *G* as reliable measurements. A modified *H* and *LE* can be calculated as:

$$LE^\* = \frac{(\mathbb{R}\_{\mathbb{N}} - G)}{(\mathbb{B} + 1)}\tag{7}$$

$$H^\star = \frac{(R\_{\mathbb{R}} - G)}{\left(\frac{1}{B} + 1\right)}\tag{8}$$

where *LE\** and *H\** denotes the closure adjusted latent and sensible heat flux, respectively.

#### *2.4. Contextual Spatial Domain*

The representative *TSEB2T* modeling grid size for the vineyard blocks was taken at 3.6 m, which corresponds to encompassing 6 × 6 grid or 36 *sUAS* thermal pixels having a resolution of 0.6 m. At this grid size, the inputs to *TSEB2T* incorporate the thermal-IR and optical bands of a vine row and adjusted interrows having a length scale of 3.35 m. Larger spatial domain scales were considered in this study, including 7.2 m, 14.4 m, and 30 m, to investigate the influence of domain size on the *TSEB2T* estimates. These selected values correspond to multiple vine rows spacing of 7.2 m (two rows), 14.4 m (four rows), and 30 m (Landsat scale—nine rows).

#### *2.5. TSEB2T Inputs*

The *TSEB2T* model developed by Nieto et al. [2] and implemented in Python language and is available at https://github.com/hectornieto/pyTSEB.

#### 2.5.1. Leaf Area Index (*LAI*)

To assess the spatial heterogeneity of *LAI*, an approach was developed in this study to calculate *LAI* using a genetic programming (*GP*) model using the Eureqa software. The *GP* model associated *sUAS* imagery and *LAI* ground measurements collected with an indirect method using (*LAI-2200*, *LI-COR*, Lincoln, Nebraska) plant canopy analyzer measurements at several locations within the northern and southern vineyards with additional validation using destructive *LAI* sampling at several locations [46]. Before performing the *GP* model calculations, imagery features were classified into two categories, vine and interrow, and then statistical calculations were separately carried out for the optical properties of each category. The main optical reflectance used in this analysis comprise the original bands (red (*R*), green (*G*), blue (*B*), and near-infrared (*NIR*)), along with two conventional *VIs* (*NDVI* and *NIR*/*R*). Statistical computations were performed using the fine-resolution data inside the spatial domain scales (3.6 m, 7.2 m, 14.4 m, and 30 m), which included the maximum, minimum, area, mean, standard deviation, and sum. The *GP* model integrates all of these corresponding statistics to construct a relationship to *LAI* observations.

### 2.5.2. Canopy Height (*hc*)

Spatial data from the digital terrain model (*DTM*) [67] and digital surface model (*DSM*) were aggregated into multiple spatial scales by employing a simple averaging method; then, *h<sup>c</sup>* was calculated using the expression: *h<sup>c</sup>* = *DSM-DTM*. For example, in the case of a 7.2-m grid, the average values of *DSM* and *DTM*, *DSM*(7.2), and *DTM*(7.2), respectively, were computed inside the grid window, then the height of the canopy was computed as: *hc*(7.2) = *DSM*(7.2) − *DTM*(7.2) .
