A Non-Reference Temperature Histogram Method for Determining Tc from Ground-Based Thermal Imagery of Orchard Tree Canopies

Obtaining average canopy temperature (Tc) by thresholding canopy pixels from on-ground thermal imagery has historically been undertaken using ‘wet’ and ‘dry’ reference surfaces in the field (reference temperature thresholding). However, this method is extremely time inefficient and can suffer inaccuracies if the surfaces are non-standardised or unable to stabilise with the environment. The research presented in this paper evaluates non-reference techniques to obtain average canopy temperature (Tc) from thermal imagery of avocado trees, both for the shaded side and sunlit side, without the need of reference temperature values. A sample of 510 thermal images (from 130 avocado trees) were acquired with a FLIR B250 handheld thermal imaging camera. Two methods based on temperature histograms were evaluated for removing non-canopy-related pixel information from the analysis, enabling Tc to be determined. These approaches included: 1) Histogram gradient thresholding based on temperature intensity changes (HG); and 2) histogram thresholding at one or more standard deviation (SD) above and below the mean. The HG method was found to be more accurate (R2 > 0.95) than the SD method in defining canopy pixels and calculating Tc from each thermal image (shaded and sunlit) when compared to the standard reference temperature thresholding method. The results from this study present an alternative non-reference method for determining Tc from ground-based thermal imagery without the need of calibration surfaces. As such, it offers a more efficient and computationally autonomous method that will ultimately support the greater adoption of non-invasive thermal technologies within a precision agricultural system.


Introduction
In agriculture, the average canopy temperature (T c ) extracted from digital thermal imagery (TI) is increasingly being used to non-destructively quantify plant water status and therefore help in the more efficient irrigation of crops [1][2][3][4][5]. In order to achieve this, thermal images have to be pre-processed to exclude non-leaf pixels (tree stems, shiny leaves, background sky, and other material) depending on the angle of view [6,7].
Currently, specialised image analysis software that often accompanies the thermal cameras, e.g., FLIR QuickReport®, Reporter Pro®, and FLIR ThermaCAM™ Researcher software (FLIR Systems, USA) [7,8], can be used to undertake a limited amount of thermal image pre-processing. However, Figure 1. Flow diagram of methodologies used in the study; HG-histogram gradient thresholding, SD-standard deviation envelope thresholding, and RT-reference temperature thresholding.
A camera-to-tree distance of 4 m was maintained for all sampled trees to produce comparable 320 × 240 pixel resolution JPEG images. A camera elevation angle of 45 0 above the horizontal was adopted to ensure only canopy features with clear sky (no sun disk) were included in the background of each image. Due to the orientation of tree rows (north-south direction) there was a high possibility of including sun disk in the camera view, which could distort the canopy image. This view angle also ensured that tree canopies from adjacent rows behind the selected target trees were not included [15]. All images were acquired during mid-day from 1100 h to 1400 h (AEST) with minimal wind conditions and clear skies, conditions that coincide with the maximum evapotranspiration demand of the canopies [7]. Relative humidity and air temperature at mid-canopy level was measured using a hygrometer throughout the image survey (Analog P2501, GH Zeal Ltd., England).
A camera-to-tree distance of 4 m was maintained for all sampled trees to produce comparable 320 × 240 pixel resolution JPEG images. A camera elevation angle of 45 • above the horizontal was adopted to ensure only canopy features with clear sky (no sun disk) were included in the background of each image. Due to the orientation of tree rows (north-south direction) there was a high possibility of including sun disk in the camera view, which could distort the canopy image. This view angle also ensured that tree canopies from adjacent rows behind the selected target trees were not included [15]. All images were acquired during mid-day from 1100 h to 1400 h (AEST) with minimal wind conditions and clear skies, conditions that coincide with the maximum evapotranspiration demand of the canopies [7]. Relative humidity and air temperature at mid-canopy level was measured using a hygrometer throughout the image survey (Analog P2501, GH Zeal Ltd., England).

Reference Surfaces
To calibrate and then validate the non-reference methods, thermal images of 'dry' and 'wet' reference surfaces were acquired during the image survey. These reference surfaces were prepared according to Jones et al. [6] using four mature leaves detached from the canopy. Two leaves serving as 'dry' references had a layer of liquid petroleum jelly (Vaseline) applied to both the abaxial and adaxial sides 1 minute before image acquisition. The two 'wet' reference leaves had a solution of water and detergent 0.01 % (v/v) applied. These references were then clipped onto a flat panel attached to a steel telescopic pole (>4 m) and positioned beside the tree canopy so as to ensure that they were included in the field of view (FOV) of each respective thermal image (Figure 2a). Precise orientation of tree rows enabled the capture of references on both shaded and sunlit sides of the trees. Thermal images of reference surfaces were only obtained every 3-4 trees, as opposed to every individual tree. The 'wet' and 'dry' values for those trees that did not have reference surfaces in the FOV were inferred from interpolation of those nearest tree values from either side of them, as outlined below. To calibrate and then validate the non-reference methods, thermal images of 'dry' and 'wet' reference surfaces were acquired during the image survey. These reference surfaces were prepared according to Jones et al. [6] using four mature leaves detached from the canopy. Two leaves serving as 'dry' references had a layer of liquid petroleum jelly (Vaseline) applied to both the abaxial and adaxial sides 1 minute before image acquisition. The two 'wet' reference leaves had a solution of water and detergent 0.01 % (v/v) applied. These references were then clipped onto a flat panel attached to a steel telescopic pole (>4 m) and positioned beside the tree canopy so as to ensure that they were included in the field of view (FOV) of each respective thermal image (Figure 2a). Precise orientation of tree rows enabled the capture of references on both shaded and sunlit sides of the trees. Thermal images of reference surfaces were only obtained every 3-4 trees, as opposed to every individual tree. The 'wet' and 'dry' values for those trees that did not have reference surfaces in the FOV were inferred from interpolation of those nearest tree values from either side of them, as outlined below.

Thermal Image Analysis
Obtained thermal images of canopies and reference surfaces were analysed separately to derive surface temperatures (degrees Celsius). Canopy thermal images were bulk processed by a customcoded script [16] to output a temperature data file and construct a temperature histogram for each image. Thermal images of reference surfaces were processed manually by using commercially available software FLIR ThermaCAM ™ researcher (FLIR Systems, USA) to obtain temperatures TwetRef and TdryRef for each side of each tree (namely sunlit and shaded) [7].

Canopy Images: Generating Temperature Histograms
In order to generate a temperature histogram for each image, the 510 FLIR JPEG images (pixels 320 × 240) of the 130 tree canopies (Figure 3a) were initially saved as 16 bit PNG images before the

Thermal Image Analysis
Obtained thermal images of canopies and reference surfaces were analysed separately to derive surface temperatures (degrees Celsius). Canopy thermal images were bulk processed by a custom-coded script [16] to output a temperature data file and construct a temperature histogram for each image. Thermal images of reference surfaces were processed manually by using commercially available software FLIR ThermaCAM ™ researcher (FLIR Systems, USA) to obtain temperatures T wetRef and T dryRef for each side of each tree (namely sunlit and shaded) [7].

Canopy Images: Generating Temperature Histograms
In order to generate a temperature histogram for each image, the 510 FLIR JPEG images (pixels 320 × 240) of the 130 tree canopies (Figure 3a) were initially saved as 16 bit PNG images before the data values were extracted as image digital numbers (DN) using 'Exiftool' in R (RStudio 1.1, Inc., Austria) [17]. The DN images were then rectified using the ImageMagick®7.0 (ImageMagick Studio, LLC) function in the command line from R, and then converted to temperature (degrees Celsius), assuming a surface emissivity of 0.98 for natural vegetation [18] and applying the measured object distance, relative humidity, and ambient temperature. Images were processed in R [16] using the 'EBImage' package [19] with a 2D matrix (200 × 150) of temperature values in .CSV file format generated for each image to support further processing (Figure 3b).
A custom code [16] written in MATLAB read each .CSV file and constructed histograms by re-organising the temperature pixel values from each file into equal-interval class bins (one degree Celsius) and calculating pixel number % (number of temperature pixels as a percentage of total image pixels) for each class. The 1 • C interval of class bins was considered as optimal to avoid over fitting the histogram shape and smoothing of important features [20]. Histogram values were normalised to 0% to 100% and the x-axis range was standardised to −30 to 70 • C for all the images. The resultant temperature histograms were typically unimodal ( Figure 3c). The peak (mode) and main body of the histogram was assumed to be associated with the temperature profile of the tree canopy pixels, whilst lower and upper 'tails' of the histogram were assumed to correspond to the sky and non-leaf (namely non-transpiring) material pixels, respectively (Figure 3a,c).

Reference Surfaces: Deriving Temperatures T wetRef and T dryRef
The average surface temperatures T wetRef and T dryRef , from the 'wet' and 'dry' references, were extracted using FLIR ThermaCAM ™ researcher Pro 2.10 software (FLIR Systems, USA). Polygons were manually drawn within the extents of the reference leaves, and emissivity was set to 0.98 as the surfaces were of natural leaves [18], an approach consistent with that reported by Fuentes et al. [7] and Jones et al. [6] (Figure 2b). As these images were acquired in sequence, a time-based interpolation between these images was carried out to derive corresponding T wetRef and T dryRef for those trees that were imaged without the referenced panels included in the FOV. The time-based interpolation was used to account for any changes in evapotranspiration demand (for example, as effected by ambient temperature, relative humidity, and wind speed) that may have occurred to the trees during the image acquisition windows. To achieve this, a polynomial curve was fitted (poly 1 or 2) to existing referenced temperature values with time, using the 'polyfit()' function in MATLAB ® 2017a (MathWorks Inc., USA) and the T wetRef and T dryRef values extracted from the curve corresponding to the image acquisition time of those non-referenced trees. Remote Sens. 2019, 11, x FOR PEER REVIEW 6 of 15

Canopy Pixels Thresholding and Obtaining Tc
The value of Tc was obtained from the post-processed images using three methods, one of them being the conventionally (and widely) used standard reference method (using physical reference panels) [7]. Each of the methods was adopted using custom code [16] written in MATLAB ® 2017a (MathWorks Inc., USA). These methods were as follows: 1. Reference temperature (RT) thresholding: Using the derived temperature values TwetRef and TdryRef from the reference surfaces, all pixel values with temperatures lower than the TwetRef or

Canopy Pixels Thresholding and Obtaining T c
The value of T c was obtained from the post-processed images using three methods, one of them being the conventionally (and widely) used standard reference method (using physical reference panels) [7]. Each of the methods was adopted using custom code [16] written in MATLAB ® 2017a (MathWorks Inc., USA). These methods were as follows:

1.
Reference temperature (RT) thresholding: Using the derived temperature values T wetRef and T dryRef from the reference surfaces, all pixel values with temperatures lower than the T wetRef or exceeding T dryRef were removed [6,7], and from the remaining pixel values, the population average canopy temperature (T c ) was calculated for each tree. The T c value for each tree (sunlit and shaded sides) was retained as the calibration/validation standard for the other two processing methods.

2.
Histogram gradient (HG) thresholding: This non-reference method is based on the identification of canopy upper and lower points (roughly corresponding to the T dry and T wet values, respectively) based on a pre-defined local gradient of the temperature histogram (Equation 1 and Figure 4a). Once assigned, the local gradient values, hitherto defined as 'ratio pixel change-RPC', were used to determine the corresponding temperature limits (RPC lower → T minRPC ; RPC upper → T maxRPC ) from which a value of T c could be derived from the population average of the pixel temperatures between the two end values (Figure 4b). An iterative process was used to symmetrically and incrementally change the predefined gradient value from 0 %/ • C (corresponding to the flat portions of the histogram tails) to the maximum gradient value (~50 %/ • C). At each pair of RPC values, the value of T c (T c-predicted ) was calculated as the population average and compared to the actual T c value (T c-actual ) derived from the standard RT method above.
Ratio pixel change (RPC) = ∆ pixel number % Class bin interval (1) 3. Standard deviation (SD) envelope thresholding: This is a non-reference method based on the standard deviation envelope filtering method of Sepúlveda-Reyes et al. [5]. In each temperature data file, the canopy T min (T wet ) and T max (T dry ) were assumed to be a numeric multiplier (m) of the standard deviation (SD) above and below the unfiltered mean temperature (T nf ), respectively.
The multiplier values ranged from m = ±0.0 (T nf ) through to ±2.0 (corresponding to the maximum and minimum values of the temperature file). Again, T c-predicted was calculated as the population average value between T min and T max for each number m and compared to the standard RT method above. and Figure 4a). Once assigned, the local gradient values, hitherto defined as 'ratio pixel change-RPC', were used to determine the corresponding temperature limits (RPClower → TminRPC; RPCupper → TmaxRPC) from which a value of Tc could be derived from the population average of the pixel temperatures between the two end values (Figure 4b). An iterative process was used to symmetrically and incrementally change the predefined gradient value from 0 %/ 0 C (corresponding to the flat portions of the histogram tails) to the maximum gradient value (~50 %/ 0 C). At each pair of RPC values, the value of Tc (Tc-predicted) was calculated as the population average and compared to the actual Tc value (Tc-actual) derived from the standard RT method above.
Ratio pixel change (RPC)= Δ pixel number % Class bin interval 3. Standard deviation (SD) envelope thresholding: This is a non-reference method based on the standard deviation envelope filtering method of Sepúlveda-Reyes et al. [5]. In each temperature data file, the canopy Tmin (Twet) and Tmax (Tdry) were assumed to be a numeric multiplier (m) of the standard deviation (SD) above and below the unfiltered mean temperature (Tnf), respectively.
The multiplier values ranged from m = ±0.0 (Tnf) through to ±2.0 (corresponding to the maximum and minimum values of the temperature file). Again, Tc-predicted was calculated as the population average value between Tmin and Tmax for each number m and compared to the standard RT method above.

Statistical Analysis
Each of the two non-reference methods (HG and SD) required iterations to determine the optimum value of the conditioning parameter, namely RPC or m, respectively. In order to determine these optimum values, 357 of the 510 images were used to 'calibrate' the methods insofar as producing T c-predicted as a function of RPC or m for comparison with the reference-derived T c-actual values of the same trees. The optimum values of RPC or m were then identified on the basis of the minimum root mean square error (RMSE). Furthermore, to select the best method of predicting T c , T c-predicted was derived again by each method (HG and SD) at their optimum RPC and m, respectively, and a linear regression analysis was performed with the corresponding T c-actual values (using 'fitlm()' function in Curve Fitting Toolbox ® in MATLAB), from which R 2 and root mean square error (RMSE) were derived. The best method with strong correlation (R 2 ) and minimum RMSE was then further validated using the remaining 153 images.

Histogram Gradient Determination of T c
The result of iterating the RPC values from 0 to 50.0 %/ • C on the derived value of T c (T c-predicted ) is depicted in Figure 5a,b for shaded and sunlit sides of the tree canopy, respectively (only RPC range 0 to 5.0 %/ • C is shown for clear illustration). Here, the RMSE and R 2 values are plotted as a function of the tested RPC values. It is noteworthy that in both Figure 5a The lower-valued minima of each graph correspond to temperature bin numbers closer to the extremes of the histograms, namely the blue sky (lower end) and non-photosynthetic material/sun-glint (higher end). Given the likelihood that these pixels are non-canopy related (allowing for sun-glint from the shiny leaf surfaces to be included here), then the upper of the two minima should be used in each case. Consequently, for the shaded-side images, the lowest RMSE (T c ) (±0.6 • C) and highest R 2 (0.83) were reported when T c was obtained at RPC 1.0 %/ • C (Figure 5a). For the sunlit-side images, the lowest RMSE (T c ) (±0.7 • C) and highest R 2 (0.83) occurred for T c obtained at RPC 0.5 %/ • C (Figure 5b). Remote Sens. 2019, 11, x FOR PEER REVIEW 9 of 15 An examination of the scatterplots (Figure 6 and 7) of Tc-actual versus Tc-predicted for each of the optimum values of RPC (shaded: 1.0 %/ 0 C; sunlit: 0.5 %/ 0 C) reveals a strong linearity (R 2 = 0.83) although there are a small number of prominent outliers (marked in red).
Outliers were further evaluated by visually observing the threshold limits obtained by each method in FLIR ThermaCAM ™ researcher Pro 2.10 software (FLIR Systems, USA) ( Figure 8). It was clear that the RT method had performed the incorrect thresholding by excluding a large proportion of the canopy pixels (white area) (Figure 8a), whereby the HG method canopy pixels were correctly thresholded (Figure 8c). This error reported by the RT method might be due to sudden environmental changes (e.g., wind, broken sunlight) which may result in alteration of actual Tc. The removal of these outliers increased the R 2 to 0.96 and 0.97 in both shaded and sunlit sides, respectively (Figure 6 and 7). An examination of the scatterplots (Figures 6 and 7) of T c-actual versus T c-predicted for each of the optimum values of RPC (shaded: 1.0 %/ • C; sunlit: 0.5 %/ • C) reveals a strong linearity (R 2 = 0.83) although there are a small number of prominent outliers (marked in red).
Outliers were further evaluated by visually observing the threshold limits obtained by each method in FLIR ThermaCAM ™ researcher Pro 2.10 software (FLIR Systems, USA) ( Figure 8). It was clear that the RT method had performed the incorrect thresholding by excluding a large proportion of the canopy pixels (white area) (Figure 8a), whereby the HG method canopy pixels were correctly thresholded (Figure 8c). This error reported by the RT method might be due to sudden environmental changes (e.g., wind, broken sunlight) which may result in alteration of actual T c . The removal of these outliers increased the R 2 to 0.96 and 0.97 in both shaded and sunlit sides, respectively (Figures 6 and 7).

Standard Deviation-Based Threshold Determination of T c
From the result of iterating the m values from ±0.0 to ±2.0 on the derived value of T c (T c-predicted ), at m = ±1.5 a minimum RMSE was recorded for both shaded-and sunlit-side images (results not shown). When examining the scatterplots (Figure 9a,b) of T c-actual versus T c-predicted for the optimum value of m = 1.5, poor correlations were revealed for both shaded and sunlit sides (shaded: R 2 = 0.02; sunlit: R 2 = 0.07) with a large number of outliers (marked in red).

Standard Deviation-Based Threshold Determination of Tc
From the result of iterating the m values from ±0.0 to ±2.0 on the derived value of Tc (Tc-predicted), at m = ±1.5 a minimum RMSE was recorded for both shaded-and sunlit-side images (results not shown). When examining the scatterplots (Figure 9a and 9b) of Tc-actual versus Tc-predicted for the optimum value of m = 1.5, poor correlations were revealed for both shaded and sunlit sides (shaded: R 2 = 0.02; sunlit: R 2 = 0.07) with a large number of outliers (marked in red). It was found that the outliers corresponded with thermal images with partial canopy inclusions ( Figure 10). When comparing Figures 10a, 10b, and 10c, during a smaller canopy situation, thresholding was poorly performed by the SD method (exclusion of canopy: White area in Figure  10c). In the cases of large canopy situations (Figures 10d, 10e, and 10f), the SD method was able to correctly threshold the canopy pixels (Figure 10f). The removal of these outliers increased the R 2 to 0.82 in both shaded and sunlit sides, respectively (Figures 9a and 9b).  It was found that the outliers corresponded with thermal images with partial canopy inclusions ( Figure 10). When comparing Figure 10a-c, during a smaller canopy situation, thresholding was poorly performed by the SD method (exclusion of canopy: White area in Figure 10c). In the cases of large canopy situations (Figure 10d-f), the SD method was able to correctly threshold the canopy pixels (Figure 10f). The removal of these outliers increased the R 2 to 0.82 in both shaded and sunlit sides, respectively (Figure 9a,b).

Standard Deviation-Based Threshold Determination of Tc
From the result of iterating the m values from ±0.0 to ±2.0 on the derived value of Tc (Tc-predicted), at m = ±1.5 a minimum RMSE was recorded for both shaded-and sunlit-side images (results not shown). When examining the scatterplots (Figure 9a and 9b) of Tc-actual versus Tc-predicted for the optimum value of m = 1.5, poor correlations were revealed for both shaded and sunlit sides (shaded: R 2 = 0.02; sunlit: R 2 = 0.07) with a large number of outliers (marked in red). It was found that the outliers corresponded with thermal images with partial canopy inclusions ( Figure 10). When comparing Figures 10a, 10b, and 10c, during a smaller canopy situation, thresholding was poorly performed by the SD method (exclusion of canopy: White area in Figure  10c). In the cases of large canopy situations (Figures 10d, 10e, and 10f), the SD method was able to correctly threshold the canopy pixels (Figure 10f). The removal of these outliers increased the R 2 to 0.82 in both shaded and sunlit sides, respectively (Figures 9a and 9b).  A summary of the statistical results obtained for both shaded and sunlit sides by HG and SD methods are shown in Table 1. HG-histogram gradient method, RPC-ratio of pixels change, SD-standard deviation, T nf -unfiltered mean temperature The HG method was found to be far more accurate than the SD method in the prediction of T c for both the shaded (R 2 = 0.83) and sunlit (R 2 = 0.83) sides of the tree (Table 1).

Validation Results
The results from the calibration work indicated the potential of the newly proposed HG thresholding method for predicting T c . To validate this outcome, the analysis was repeated on additional canopy images. According to Figure 11a,b, strong linearity was observed for both shaded (R 2 = 0.90) and sunlit (R 2 = 0.89) side images, similarly to the calibration results. There were few outliers resulting from the error in T c_actual (marked in red), where removal of them improved the correlation to R 2 = 0.96 for both sides (Figure 11a,b). A summary of the statistical results obtained for both shaded and sunlit sides by HG and SD methods are shown in Table 1. The HG method was found to be far more accurate than the SD method in the prediction of Tc for both the shaded (R 2 = 0.83) and sunlit (R 2 = 0.83) sides of the tree (Table 1).

Validation Results
The results from the calibration work indicated the potential of the newly proposed HG thresholding method for predicting Tc. To validate this outcome, the analysis was repeated on additional canopy images. According to Figures 11a and 11b, strong linearity was observed for both shaded (R 2 = 0.90) and sunlit (R 2 = 0.89) side images, similarly to the calibration results. There were few outliers resulting from the error in Tc_actual (marked in red), where removal of them improved the correlation to R 2 = 0.96 for both sides (Figures 11a and 11b).
Therefore, the resulting 'optimum thresholds' for the HG method were found to be consistent with the threshold temperature histogram at RPC 1.0 %/ 0 C and 0.5 %/ 0 C for shaded and sunlit images, respectively. This was true even for different canopy sizes.

Discussion
In order to determine the average temperature of a tree canopy (Tc), it is necessary to remove non-canopy-related pixels, such as stem parts, sky, and leaves exhibiting specular reflection from thermal imagery. Historically, the inclusion of both TwetRef and TdryRef within each thermal image has been used to achieve this thresholding. Although effective, this method is time consuming and labour intensive, and therefore prohibitive to the collection of thermal information for many trees, i.e., commercial orchards at any one time. The histogram gradient (HG) method presented in this study offers a histogram method for thresholding canopy pixels from thermal images that produces results comparable to those from the accepted reference temperature (RT) thresholding method. The Therefore, the resulting 'optimum thresholds' for the HG method were found to be consistent with the threshold temperature histogram at RPC 1.0 %/ • C and 0.5 %/ • C for shaded and sunlit images, respectively. This was true even for different canopy sizes.

Discussion
In order to determine the average temperature of a tree canopy (T c ), it is necessary to remove non-canopy-related pixels, such as stem parts, sky, and leaves exhibiting specular reflection from thermal imagery. Historically, the inclusion of both T wetRef and T dryRef within each thermal image has been used to achieve this thresholding. Although effective, this method is time consuming and labour intensive, and therefore prohibitive to the collection of thermal information for many trees, i.e., commercial orchards at any one time. The histogram gradient (HG) method presented in this study offers a histogram method for thresholding canopy pixels from thermal images that produces results comparable to those from the accepted reference temperature (RT) thresholding method. The development of the unimodal histogram for each image and the subsequent removal of the left and right 'tails' (e.g., thresholding limits T min and T max of the canopy) successfully removed pixel information associated with non-canopy-related features. The consistency of histograms derived from each image was maintained by including the sky, without the sun disc or canopies from surrounding trees, in each image. The consistent angle of view and 4 m distance between the camera and the tree canopy also supported a standardised spatial scaling.
The FLIR B250 camera has a low IR resolution of 200 × 150 pixels, which may result in a lesser number of pixels of the target than could be obtained from a high-resolution camera. However, the proposed method should be applicable even across different camera resolutions as far as the histogram is constructed with one degree Celsius class bins and the number of pixels per class bin is presented as a percentage of the total number of pixels in the image.
Whilst Sepúlveda-Reyes et al. [5] achieved encouraging results using the SD method for thresholding aerial thermal images of tree canopies, the same result was not achieved in this study using ground-based thermal images. This may be due to the temperature histograms from ground thermal images being skewed left (Figure 3c) due to the inclusion of lower-temperature pixels from the sky. Furthermore, a Gaussian mixture method was previously suggested by Wang et al. [3] to derive T c from ground thermal images, but it required the temperature pixel distribution to be normalised, which was not required in this HG method.
The limitations of the new HG thresholding method are that it requires the maintenance of consistency during image acquisition (viewing angle, distance to tree, and clear sky), and the background temperature must be distinguishable (either low or high) from canopy temperature. As a future improvement, additional non-reference procedures based on histogram thresholding could be tested. This may include Gaussian mixture distribution extraction [3,21], colour intensity histograms [22], or using deep machine learning algorithms such as convolution neural networks [23]. Additionally, in future the performance of the new method should be evaluated over tree species other than avocado, different environmental situations such as early morning and late afternoon, cloudy conditions, and on different platforms such as mounted on a ground vehicle However, the histogram gradient threshold presents as a new and easy way of thresholding canopy pixels in ground canopy thermal imagery, enabling the monitoring of individual trees for changes in water status caused by abiotic and biotic factors. Therefore, this method offers great benefit to future agriculture processes (precision farming) that rely on large volumes of thermal remotely sensed imagery.

Conclusions
The results from this study present an alternative, histogram method for thresholding canopy area from ground-based thermal imagery without the need of temperatures from reference surfaces. By developing a unimodal histogram for each image and then removing the left and right 'tails' associated with non-canopy-related features based on histogram gradients (e.g., T wet and T dry boundaries) predictions of average canopy temperature (T c ) were achieved. These predictions of T c were found to be highly comparable to T c produced from the commonly used method (reference temperature thresholding), which requires the inclusion of physical 'wet' and 'dry' references within each image. This novel 'histogram gradient thresholding' method offers a more time-efficient and labour-saving method, and is a potentially automated method of calculating canopy temperature from thermal imagery, allowing large volumes of images to be processed in a short time. This result was not previously possible with either the reference method or the manual digitisation of canopy specific pixels within each image. With the ever increasing interest in remote sensing technologies for monitoring tree health, this outcome offers significant benefit for the rapid and standardised analysis of thermal imagery.