Monitoring LAI, Chlorophylls, and Carotenoids Content of a Woodland Savanna Using Hyperspectral Imagery and 3D Radiative Transfer Modeling

: Leaf pigment contents, such as chlorophylls a and b content (C ab ) or carotenoid content (Car), and the leaf area index (LAI) are recognized indicators of plants’ and forests’ health status that can be estimated through hyperspectral imagery. Their measurement on a seasonal and yearly basis is critical to monitor plant response and adaptation to stress, such as droughts. While extensively done over dense canopies, estimation of these variables over tree-grass ecosystems with very low overstory LAI (mean site LAI < 1 m 2 /m 2 ), such as woodland savannas, is lacking. We investigated the use of look-up table (LUT)-based inversion of a radiative transfer model to retrieve LAI and leaf C ab and Car from AVIRIS images at an 18 m spatial resolution at multiple dates over a broadleaved woodland savanna during the California drought. We compared the performances of different cost functions in the inversion step. We demonstrated the spatial consistency of our LAI, C ab , and Car estimations using validation data from low and high canopy cover parts of the site, and their temporal consistency by qualitatively confronting their variations over two years with those that would be expected. We concluded that LUT-based inversions of medium-resolution hyperspectral images, achieved with a simple geometric representation of the canopy within a 3D radiative transfer model (RTM), are a valid means of monitoring woodland savannas and more generally sparse forests, although for maximum applicability, the inversion cost functions should be selected using validation data from multiple dates. Validation revealed that for monitoring use: the normalized difference vegetation index (NDVI) outperformed other indices for LAI estimations (root mean square error (RMSE) = 0.22 m 2 /m 2 , R 2 = 0.81); the band ratio ρ 0.750 µ m ρ 0.550 µ m retrieved C ab more accurately than other chlorophylls indices (RMSE = 5.94 µ g/cm 2 , R 2 = 0.75); RMSE over the 0.5–0.55 µ m interval showed encouraging results for Car estimations.


Introduction
The Mediterranean-like ecoregions are expected to face the greatest proportional change in their biodiversity by 2100 because of climate change [1]. These ecoregions (Mediterranean Basin, California, Chile, South Africa and Western Australia) are characterized by hot and dry summers, cool and wet winters, and a large annual rainfall variability. They present a high plant and animal diversity spread within numerous vegetation types [2]. Mediterranean forests in particular range from dense forests to dehesa pastures, with trees scattered over a grass field. They notably display a high proportion of broadleaved species compared to forests found in northern parts of Europe [2]. While plants have adapted to this climate, the Mediterranean biome is still one of the most imperilled, as increased frequency of heat and drought periods as well as growing anthropogenic pressure are threatening the present equilibrium [3][4][5]. To identify key threats and guide forest conservation efforts, it is crucial to understand disturbance regimes, species adaptability, and ecological processes [6]: An improved monitoring of the ecoregions is therefore necessary. Essential biodiversity variables (EBVs) have been recently conceptualized by the Biodiversity Observation Network of the Group on Earth Observations (GEOBON) to unify the global monitoring efforts [7,8]. Foliar pigment contents (such as chlorophylls a + b content C ab , carotenoids content Car) and leaf area index (LAI) have been identified as traits associated with the ecosystem structure and its relationship to biodiversity. Indeed, leaf pigments are involved in the photosynthetic process, leaf protection, and also contribute to the ecosystem's biogeochemical and nutrient cycles [9][10][11], while LAI controls many biological and physical processes such as gross primary production (GPP), rainfall interception, and heat fluxes [12]. Estimating these variables should therefore be a priority to guide conservation policies, as changes in the temporal patterns of LAI and pigment variations could be due to vegetation stress.
Field sampling and laboratory measurements of these traits are not acceptable methods to monitor vegetation, as they are time-consuming, spatially limited, and cannot capture inter-and intraspecies variability. Remote sensing, using either multispectral or hyperspectral sensors, has been extensively used to map LAI and leaf pigments over large swathes [13][14][15][16] and allows for recurrent revisits of the study sites. Spectral images are obtained by sensors measuring forest canopy reflectance over a large spectrum, so that slight variations can be detected in the forest canopy's reflected radiation. However, these variations are not exclusively due to the vegetation. Indeed, illumination, local atmospheric conditions, and ground reflectance also have an influence on the measured reflectance [17][18][19]. Ground reflectance contribution to the measured spectrum can be most problematic for open canopies or sparse forests (which have very low LAI) such as woodland savannas, which are one of the major plant communities in the Mediterranean biome. Indeed, while at high spatial resolutions vegetation pixel could be clearly identified, reducing ground influence, pixels from images with medium to low spatial resolution will only have a limited canopy contribution to the overall reflectance.
Remote sensing pigment and LAI estimation methods can belong to three main families: empirical-statistical approaches, which fit a model to a set of validation data acquired in-field; physical approaches, which rely on radiative transfer model (RTM) inversions; and hybrid methods, which are still in their infancy and combine the properties of physically-based methods with the flexibility of the empirical-statistical approaches [20]. While empirical-statistical approaches have the disadvantage of difficult transferability to other sensors, sites, or dates, as the model may be overfit, RTM inversions' principal drawback is their ill-posedness: An RTM simulates top-of-canopy reflectance using a modeled scene and physical principles to compute the radiative budget, and various sets of inputs may lead to similar outputs as the problem faces too many unknown variables. However, various options have been successfully explored to mitigate these issues and obtain reliable solutions, such as using prior information to lower the variability of the estimations' distribution or introducing spatial constraints during the inversion process (see an overview of possible options in Baret and Buis [21]). In the present study, a look-up table (LUT)-based RTM inversion was performed using the DART 3D canopy model [22], based on the work done in Adeline et al. [23].
The aim of this study was to find an inversion strategy suited for a heterogeneous tree-grass ecosystem's health (e.g., resiliency to change, continuity in GPP, and nutrient capital [24]) monitoring, i.e., capable of estimating LAI and leaf C ab and Car over time using multitemporal hyperspectral airborne acquisition with medium ground sampling distance.
Canopy LAI and biochemistry estimations using hyperspectral imagery have mostly been done over canopies with high LAI (for instance, see Banskota et al. [16], le Maire et al. [25], Ali et al. [26], Darvishzadeh et al. [27], Malenovský et al. [28]). While various estimations have also specifically been done over open-canopy ecosystems (e.g., the works of Zarco-Tejada et al. [29], Hernández-Clemente et al. [30], Zarco-Tejada et al. [31]), research concerning acceptable modeling methods within RTM is still ongoing [32,33]. The validity of a simplified representation of trees (simple forest representation, SFR) within the DART scene (such as the one done by Gascon et al. [14]) when working with medium-resolution hyperspectral images of very sparse forests was tested in this study, based on the work presented in Widlowski et al. [32]. As acceptable representation of coniferous trees within RTM is still a difficult topic [34]; this study only focused on broadleaved trees as a preliminary test case to study Mediterranean forests. Validation data were available at multiple dates, allowing for a more robust calibration of the inversion method. The validity of the results presented in this paper was assessed using: (i) comparison between estimations and field measurements when and where available, and (ii) comparison between the site's expected and estimated LAI and foliar pigment contents variations over two years.

Study Site
The study site is an oak woodland savanna located in the lower foothills of the Sierra Nevada Mountains (Tonzi Ranch, latitude: 38.431 • N; longitude: 120.966 • W; altitude: 177 m; average slope: 1.5 • ). The climate is Mediterranean, with wet winters and hot, dry summers without rainfall. The overstory is mainly composed of blue oaks (Quercus Douglasii-QUDO), with a sparse (<10%) population of grey pines (Pinus Sabiniana-PISA). Blue oak trees start to become active in April and have shed most of their leaves by November. Spatial distribution and structural composition are typical of heterogeneous canopies, with a mean canopy cover (CC) of 47% and a mean LAI of 0.8 m 2 /m 2 . The understory is composed of cool season annual C3 grass species: It is active from December to May (wet period) and dry during the summer and autumn. Both oak trees and grasses can therefore be active in April and May. The soil is an Auburn very rocky silt loam (Lithic haploxerepts). More detailed site information can be found in previous studies [35][36][37]. An aerial picture of the site and its location are given in Figure 1.

Leaf Area Index
Field data were collected coincident with NASA Hyperspectral Infrared Imager (HyspIRI) Mission Study Airborne Campaigns, which took place in fall 2013 and summers 2014 and 2016 (https://hyspiri. jpl.nasa.gov/). Several digital hemispherical photographies (DHP) were collected over multiple 20 m × 20 m plots across the study site covering all vegetation cover fraction and species composition. The LAI of each plot is therefore considered as an LAI measurement. Information concerning the number of plots for each date is given in Table 1.
From each 20 m × 20 m plot, three DHP were taken from a Nikon Coolpix 4300 camera post sunset when no direct sunlight is visible. DHP were processed using Hemiview software [38]. Hemiview calculates PAI from the fisheye images estimating image gap fraction using eight azimuth and seven zenith angles. The threshold was estimated for each image based on a visual interpretation. Hemiview calculates an effective plant area index (PAI) that has to be corrected by clumping factors. However, some studies have considered that in forests, effective PAI is almost identical to true LAI, as the underestimation of LAI due to clumping effects could be compensated by the overestimation due to woody structures [39]. In the rest of the paper, effective PAI estimated through Hemiview is therefore considered equal to true LAI and simply referred to as LAI.

Leaf Biochemistry
A set of leaves from five QUDO individuals were collected from the upper, sunlit portion of the canopy. Sampling started within an hour of the timing of the overflight. Leaf samples were collected from open grown trees that were in full sunlight, as high into the canopy as possible to reach, and from branches on the east and west side of the tree. Leaves were placed in foil packets so that light would not reach them and stored in a container on blue ice and transferred to a lab refrigerator at the field station until lab measurements could be made, no later than 48 hours after data collection. Leaf samples were frozen in liquid nitrogen for a short time until lyophilized.
The analysis extraction was done in 90% acetone (14 mL) for 48 hours, which provided enough solvent extract to make 3.5 mL samples for each individual. Lambda 25 UV/Vis spectrophotometer runs were done with standards with a concentration of 23 mg/mL for chlorophylls and 8 mg/mL for carotenoids. Absorbance A of the solution was measured at 0.470, 0.662, 0.645, and 0.710 µm. Chlorophylls and carotenoids concentrations of the solution (c a , c b , c (x+c) ) were computed using Equations (1)-(3) [40,41] with a dilution factor equal to 3.5, and further processed to get chlorophylls a + b and carotenoid concentration on the leaf area (µm/cm 2 ). Information concerning the final number of validation points is given in Table 1.

Trunk Reflectances
Trunk reflectances of QUDO trees were collected using an Analytical Spectral Devices, ASD Inc., Boulder, CO, USA contact probe. The spectroradiometer was calibrated using a spectralon panel as a white reference before every acquisition. Trunk reflectances were obtained over the 0.350-2.500 µm spectral range. AVIRIS-C hyperspectral data are processed and delivered by NASA Jet Propulsion Laboratory (JPL; http://aviris.jpl.nasa.gov). The sensor is a whiskbroom scanner with an instantaneous field of view of 1 mrad and a total field of view of 34 • . It is composed of 224 contiguous bands with a full width at half maximum (FWHM) of 0.010 µm from 0.365 to 2.5 µm [42]. Airborne acquisitions took place in spring, summer, and fall (corresponding to the 3 phenological stages) from 2013 to 2015 and every summer from 2016 to 2018. Images were acquired at nadir. Acquisitions used in this study are given in Table 2. The ER-2 aircraft flew 20 km above the ground, around noon, in order to avoid spectral directional effects. Preprocessing steps provided by NASA JPL included radiometric calibration, geometrical orthorectification, nearest neighbor spatial resampling at 18 m, and atmospherical correction performed with ATREM [43], in order to retrieve surface reflectance. Additional flights over the study site were done during the summer 2014 using the AVIRIS-Next Generation (AVIRIS-NG) sensor. AVIRIS-NG images have a 2 m spatial resolution and 432 spectral bands with a FWHM of 0.005 µm from 0.380 to 2.510 µm.

General Methodology
In the present section, a description of each step of the methodology followed in this study is given. Figure 2 presents an outline of the steps described in the following subsections.

Image Processing for Multitemporal AVIRIS-C Data
In addition to the image processing protocol from NASA JPL, co-registration was further considered to improve the spatial matching between the AVIRIS-C images over the Tonzi site. The reference image was the summer 2014 AVIRIS-NG image, as the georeferencement seemed acceptable when overlaid on QGIS Microsoft Bing Maps satellite backgrounds. Co-registration was based on a nearest-neighbour resampling using around 25 ground control points per image. This amount could vary depending whether or not some features were identifiable on the AVIRIS-C images. For each image, the co-registration root mean square error (RMSE) was less than one pixel. A temporal correction was also applied on the images, serving as an in-between flight radiometric normalization in order to compare the images on the same spectral baseline. An empirical line correction method [44][45][46] was done using 9 invariant targets (3 black targets from the waters of Tahoe lake; 3 grey targets from flat sand/gravel areas; 3 white from different flat roofs with white painting) within the AVIRIS-C atmospherically corrected reflectance images. The reference image for the temporal correction was the spring 2013 AVIRIS-C image. The calibration targets were all located outside of the Tonzi site and spectral comparison over the various AVIRIS-C images showed that they were invariant in time.

Masking Using Canopy Cover and Species Composition Maps Derived From the AVIRIS-NG Image
The objective was to create a monospecific mask of QUDO. The final monospecific mask was obtained by (i) eliminating the areas with too many PISA in the images, (ii) eliminating pixels with too low tree presence (CC < 10%) and (iii) removing the pixels containing water, roads, and manmade infrastructures. The AVIRIS-NG image at 2 m spatial resolution was processed to derive a 18 m CC map as well as a species composition map. Figure 3 gives an overview of the processing steps done over the AVIRIS-NG image.
To build a canopy cover map, a vegetation mask was built based on a simple histogram thresholding with the application of a vegetation index. Version 2 of the modified chlorophyll absorption ratio index (MCARI2, [47]) was chosen given that the two main components in the image can be easily separated in summer: a healthy green tree canopy (chlorophylls peak) and dry dead grass (no chlorophyll). Finally, a spatial resampling was applied using the pixel aggregate method to come back to the 18 m spatial resolution of AVIRIS. The statistics provided from the derived canopy cover map show that the mean canopy cover is 53% over the area of interest with a concave bell curve from 0% to 90%. In addition, a peak also appears around 100% in the distribution histogram (see Figure 3). To build the species classification map, the supervised support vector machine (SVM) classification method was used through the program LibSVM [48]. One class was created for each tree species (QUDO and PISA). A total of 506 and 557 pixels were respectively selected for each class as training samples based on fieldwork data, tree structural information derived from LiDAR data collected in 2009 (PISA sample height assumed to be higher than 20 m, [49]) and visual inspection. A total of 339 spectral bands were selected to obtain the same spectral bands as AVIRIS-C. The SVM was run with the radial basis function kernel (C: 6.5; γ: 0.0055; accuracy: 93%; cross-validation correlation coefficient: 98.4). For verification purposes, the resulting classification map was also visually compared to a true-color image of the site, as the tree species have clearly distinct colors. Finally, an 18 m species composition map was created by applying the same spatial resampling as for the canopy cover map.
Pixels where QUDO crowns represented less than 80% of the total tree crowns were rejected as the PISA presence would be too high, as well as pixels with a CC below 10%. Finally, pixels containing water, roads, and manmade infrastructures were manually masked over each image.

DART and PROSPECT Radiative Transfer Models Parametrization
DART version 5.7.3v1078 was used to simulate canopy reflectances and build the LUT. DART is a radiative transfer model able to simulate light interactions and multiple scattering effects within the modeling of a 3D scene, including the topography and the atmosphere. The 3D scene is a made of voxels whose size can be set by the user to simulate various scene elements. These elements can be represented through planar elements (triangles, allowing for precise modeling) and/or turbid voxels (statistical representation of fluids and vegetation). In the case of vegetation, these turbid voxels properties are determined by the canopy leaf area and its spatial distribution, the leaf angle distribution (LAD), and leaf optical properties. A more precise description of the DART model can be found in Gastellu-Etchegorry et al. [22] and Gastellu-Etchegorry et al. [50]. DART makes it possible to define trees with simplified representations such as turbid crowns of conic, spherical, or ellipsoidal shape and allows for tuning of the fraction on full voxels within the crown as well as branches and twigs modeling. Therefore, trees are defined by structural parameters such as the shape and size of their crown and the distribution of their branches and leaf cells, as well as optical parameters with the optical properties of leaf voxels and branches.
The scene modeling done in this study is based on a simplified forest representation. A repetitive pattern of four trees with ellipsoidal crowns was used [14]. Crown size and height were set to the Tonzi site average values [36]. Branches and twigs were neglected to reduce modeling complexity and introduction of errors if not correctly set. Trunk reflectance was set to the value measured in-field. Canopy cover variations were obtained by increasing the scene dimensions. The trees were placed so that the scene bidirectional reflectance factor (BRF) is the closest to a forest BRF, following the pattern given in Gastellu-Etchegorry et al. [51] and Gascon et al. [14]. Trees and scene characteristics are given in Table 3. The forest understory was represented as a lambertian surface whose optical properties were determined by pure soil pixels selected within the AVIRIS-C images in open parts of the site, at least one pixel (18 m) away from canopy pixels and several pixels wide. These open parts were visually chosen using QGIS high-resolution Microsoft Bing Maps satellite backgrounds. Leaf optical properties were simulated through the PROSPECT model [52] implemented in DART, using the input parameters described in Table 4. The leaf structure parameter N was set to 1.8. Sun and atmosphere were the only radiation sources. The sun positions were set so that they match with the AVIRIS-C images based on site location and hours of acquisition. The simulations used a rural aerosol model with a visibility of 23 km corresponding to mid-latitude summer standard atmosphere. A total of 168 bands with a spectral bandwidth of 0.010 µm were simulated between 0.36 and 2.45 µm to cover both the visible and short wave infrared (SWIR) spectral ranges.

Inversion Stategies
This study follows a LUT-based approach consisting in finding the DART model parameters minimizing a cost function comparing measured y to computedŷ reflectances. Since the inversion problem is ill-posed, multiple sets of parameters can yield similar reflectances. A common estimation method is to consider the mean values from the q best sets of parameters as the final solution [53,54]. Various cost functions were tested: root mean square error (RMSE, Equation (4)), spectral angle (SAM, Equation (5)) and vegetation index (VI) differences (D V I , Equation (6)).
RMSE and SAM were computed using variable-specific spectral intervals. The near infrared (NIR) and SWIR were used for LAI (INT LAI, 0.8-2.45 µm), while parts of the visible range were used for C ab and Car (INT CAB, 0.5-0.75 and INT CAR, 0.5-0.55 µm respectively). D V I was computed using various indices. Concerning LAI, the NDVI was selected, as it is commonly used to estimate vegetation cover, as well as MSAVI2, which is a soil-adjusted index; for C ab , Maccioni, GM_94b, gNDVI, MCARI2, and TCARI/OSAVI were considered, since they are soil-adjusted indices specifically designed for low-LAI environments, and they were successfully used in a sparse forest context [47,[55][56][57]; concerning Car, CRI and R515/R570 were chosen, as they were also successfully used in previous studies [30,57]. Table 5 summarizes the selected methods for all investigated variables as well as the number of spectral bands used.
The number q of best-matching cases to use for the inversion was selected on a per-variable basis. The objective was to find an optimal value that would minimize the RMSE between field data and estimations when performing the LUT-based inversion. Therefore, inversions were done with increasing q values, from 1 to 400 (0.9% of the LUT). If an optimal value existed across all validation dates, it would determine the number of solutions to keep for the inversion validation step (Section 2.3.5) and then the inversion of all remaining AVIRIS-C dates available.

Validation Metrics
Validation was performed by comparing LAI, chlorophylls, and carotenoids estimates with the field measurements available using RMSE, bias, the R 2 of the predicted-vs.-measured regression line and the STDB (standard deviation of the difference between estimated points and the regression line) as cost functions.
Concerning LAI, each validation point was the average of LAI values derived from 3 hemispherical pictures taken within a 20 × 20 m plot (see Section 2.2.1). Therefore, one LAI validation point does not necessarily only contain information about the LAI of the AVIRIS-C pixel it is associated with. Figure 4 illustrates the issue: Since the AVIRIS pixels and the grid used for DHP acquisitions have no reason to match, some DHP (dots) may have been acquired out of the pixel the final LAI has been associated with (red square). Direct comparison between pixel-estimated value and validation data can therefore be inappropriate, especially for spatially heterogeneous canopies. To address this, LAI validation points were compared to a spatially smoothed LAI inversion map. The spatial smoothing was done by computing the average LAI within a sliding window with an aperture size of 3 pixels. The aperture size was set to 3 as the LAI validation plot size (20 × 20 m) roughly matches the AVIRIS-C pixel size (18 × 18 m) (Figure 4).  [63] Biochemistry validation data were obtained at the leaf scale, for a single tree in each validation pixel. There is therefore an unavoidable uncertainty, as validation is done by comparing field data at the leaf scale (single tree) to plot biochemistry (canopy scale, multiple trees). An assumption is made that tree biochemistry does not vary much spatially, and C ab and Car estimations were directly extracted from the pixels associated to the acquisition positions.  2 and 3.3.3). Changes in the site's mean LAI and foliar C ab and Car were considered significant if they were greater than the STDB of their respective variables.

Comparison between Airborne and DART-Simulated Reflectances
The assessment of LUT reflectances acceptability with regard to those measured by AVIRIS-C was done in two steps. First, DART-simulated reflectances of each LUT were plotted and compared to the reflectances measured by AVIRIS-C at the associated date (see, for example, Figure 5). All LUTs consistently overestimated canopy reflectances between 0.75 and 1.1 µm , with the mean image reflectance being either barely within or below LUT reflectances near the red-edge. Image reflectances in the visible and SWIR regions were, however, within simulation ranges. Secondly, the amount of pixels that could be explained by the LUTs was evaluated. The percentages of pixels explained by the LUTs within each image are given in Table 6. LUT reflectances encompass image reflectances at various degrees depending on the spectral intervals considered (INT LAI, INT CAB, INT CAR). On average, 54% of pixels' reflectances completely fall within the LUT reflectance boundaries for INT LAI. This goes up to 67% and 97% for INT CAB and INT CAR, respectively. When considering the vegetation indices computed from both images and LUT reflectances, pixel-derived vegetation index (VI) are within LUT-derived VI boundaries from 90% to 100% of the time (excepting summer 2017, which has one occurrence at 83%), depending on the VI considered.

Influence of the Canopy Cover on Vegetation Indices
Once the LUT was built, it was also necessary to assess the suitability of each VI to serve as a criterion for LAI and biochemistry estimations when using the LUT entries. Indeed, it is clear that the green vegetation signal is quite low. The reflectance of low CC pixels is close to that of pure ground, with only a slight red-edge increase ( Figure 6). Even at high CC, the red-edge is not very steep, and the chlorophylls reflectance peak at 0.55 µm is barely visible: It could be possible that some VIs would not be able to detect some of the variations of the variables of interest.
To do this assessment, vegetation indices were computed over the LUT reflectances. Figure 7 shows the evolution of three different VIs over the variable they are associated with, on a per-CC basis (colors). Ideally, no matter the CC, the VI should present variations when their associated variable increases, i.e., behave as a bijective function.
It was found that the indices selected for LAI estimations all started saturating at 0.4 m 2 /m 2 LAI for the 10% of CC cases (Figure 7a, red line). This behavior was not found for higher CC. Since the 10% CC category could only bring confusion due to this saturation when estimating LAI, it was decided to remove these cases from the LUT. Excluding these cases could also be done based on their lack of realism, as low CC-high LAI canopies are unrealistic. Vegetation indices dedicated to C ab and Car estimation did not show such behavior (Figure 7b,c); therefore, all CC cases were considered during the inversion step for C ab and Car.

Determination of the Number q Achieving the Best Inversion Performances
For each variable, in order to assess the uncertainty associated with determining q using only validation data from a single date, q was determined both on a date-to-date basis and on a multidate basis.
For LAI, ignoring the RMSE INT criterion which performed very poorly, all criteria seemed to have reached a plateau by q = 100 at most (see Table 7 and Figure 8a). As there was no obvious optimal q which could minimize RMSE, it was decided to consider the 100 best performing solutions for LAI estimation. While for summers 2014 and 2016, SAM INT LAI, NDVI, and MSAVI2 performed similarly, for fall 2013, only the VI-based criteria achieved acceptable performances (around 0.6 m 2 /m 2 for RMSE and SAM INT LAI, around 0.2 m 2 /m 2 for NDVI and MSAVI2).
C ab results showed a similar plateauing behavior (see Table 7 and Figure 8b). Therefore, it was decided to proceed the same way as the LAI and to consider the 300 best performing solutions for C ab estimations, since the RMSE plateau was always reached at this point. The RMSE of the criteria for summer 2013 were all rather high, with only GM_94b obtaining a RMSE below 10 µg/cm 2 . The lowest RMSEs for fall 2013 and summer 2014 were obtained with gNDVI (3.86 µg/cm 2 ) and GM_94b (2.89 µg/cm 2 ), respectively. MCARI2 performed very poorly for all dates, even more so than the RMSE and SAM criteria that at least showed acceptable performances for summer 2014. For Car, an optimal q value was found when singling out the summer 2013 data and using the RMSE INT CAR criterion. However, at the other dates, the RMSE INT CAR criterion barely reached its plateauing value by q = 400 (see Table 7 and Figure 8c) Table 8 shows LAI retrieval performances when using the mean of the q = 100 best solutions. Inversions based on NDVI and MSAVI2 differences performed equivalently in terms of RMSE (less than 0.25 m 2 /m 2 ), bias, STDB, and R 2 (above 0.80). RMSE INT LAI performed the worse for all validation metrics, followed by SAM INT LAI. Figure 9 shows that the highest LAI are slightly underestimated and vice versa, and that estimation uncertainties increase with LAI. Concerning low LAI points, only the fall 2013 points (blue) can be overestimated. No clear distinction in terms of performances is observed between other seasonal measurements. The global performance shows a R 2 of 0.8 and a total RMSE of 0.22 m 2 /m 2 .

C ab and Car Retrieval Performance Comparison
For C ab at q = 300, apart from D MCARI2 , VI differences performed better than methods based on RMSE and SAM (Table 9). GM_94b is the overall best-performing VI, with the lowest RMSE, highest R 2 , and lowest STDB (5.94 µg/cm 2 , 0.75, and 4.06 µg/cm 2 , respectively), besting even soil-adjusted VI.
When compared to field measurements, most GM_94b-estimated points are very close to the first bisector, and only one point (dark orange from summer 2013) is greatly misestimated (Figure 10a). This point is in an open-canopy part of the site (Figure 1), close to the light orange one, which is well estimated. Both the open-and closed-canopies parts of the site are well estimated for summes and fall, with C ab values ranging from 20 to 53 µg/cm 2 . The global performance is promising with a R 2 of 0.75 and a total RMSE of 5.94 µg/cm 2 . For Car, at q = 400, the quality of estimation results is not as clear: RMSE can be quite low (3 µg/cm 2 or below for CRI and RMSE INT CAR), but R 2 remains low (0.32 at most for SAM INT CAR), and no linear relationship is found between estimated and measured data (Table 9). RMSE INT CAR and CRI both present low bias (≤0.7 µg/cm 2 ) and STDB (≤2.2 µg/cm 2 ).
The RMSE INT CAR method showed the best performances overall. Apart from two points at 10 µg/cm 2 that are overestimated and one point (dark orange at 12.5 µg/cm 2 ) that is severely underestimated, both the open-and closed-canopies parts of the site are well estimated (Figure 10b). The R 2 value is 0.11, and the total RMSE is 2.57 µg/cm 2 . Table 9. Inversion method performances for C ab (q = 300) and Car (q = 400) estimations over all dates available.

LAI, C ab and Car Seasonal Monitoring
In Appendix A are given estimation maps of Tonzi summer LAI and canopy C ab and Car (LAI × foliar pigment content) from 2013 to 2018 ( Figures A1-A3, in Appendix A respectively). These maps were computed using the best VI identified in Sections 3.3.2 and 3.3.3.
Focusing on years 2013 and 2015, estimated LAI and foliar C ab decreased from spring to fall (steady decrease from 0.9 m 2 /m 2 and 38 µg/cm 2 to 0.4 m 2 /m 2 and 30 µg/cm 2 , respectively). Site LAI variations from summer to fall in 2013 and 2015 are greater than the estimations' STDB ( Figure 11), while from C ar only, summer to fall variations in 2015 are greater than the STDB (Figure 12). Concerning carotenoids, estimated variations are below the STDB and are therefore not remarkable.

Limitations of the SFR for Low-LAI Sparse Forests Modeling within DART
The forest representation chosen for this study aimed at minimizing the a priori knowledge necessary to monitor the Tonzi site. Trees modeled within the DART scenes had simplistic crown shapes and few nonphotosynthetic elements, as branches were not included. This has an impact on both LAI and C ab estimations, as woody elements' influence on canopy reflectance can be important [28,56] for heterogeneous and defoliated stands such as the Tonzi site. In particular, Malenovský et al. [28] and Widlowski et al. [32] showed that inclusion of woody elements within the DART canopies affects NIR reflectances for the resolutions considered in this study, which could explain the overestimated reflectances of the LUT in this spectral region ( Figure 5).
These mismatches between computed and measured reflectances are such that some of the inversion methods tested in this study are not appropriate. Indeed, Table 6 shows that when considering the INT LAI and the INT CAB intervals, a large proportion of the AVIRIS-C images' pixels could not be explained by the LUTs. The consequence is that VI-based inversion criteria outperform those based on spectral intervals. While some VI may use reflectances in the NIR which are overestimated by DART ( Figure 5), it appears that the overestimation effect is dampened as VI are less sensitive to differences in reflectance levels, and pixel-derived VIs are almost always within the boundaries of LUT-derived VIs.
The simplistic ground representation within the DART scenes, with a flat ground having a lambertian reflectance, may also be a source of errors, especially for the spring acquisitions, as important spatial variations of ground reflectance can be expected (see Sections 2.1 and 4.4). Unfortunately, such variation cannot be taken into account by the simple forest representation implemented in this study. A possible way to tackle this, though it would increase computation time, would be to introduce ground reflectance as another LUT parameter. However, the difficulty would be to make sure that most of the possible ground reflectances are captured by the LUT, which may not be an easy task when working with medium-to low-resolution images.

Time Dependency of the Best Performing Inversion Method
The objective was to be able to monitor the study site over multiple dates with or without available field data. As field data are usually scarce, inversion methods are often calibrated using data from only a single date. While this guarantees that the optimal method for this date has been found, the validity of extrapolations to images acquired at different times is not a given. Indeed, in Section 3.3's Table 7, both GM_94b and gNDVI indices could be identified as optimal depending on the date. However, when considering the complete dataset, which includes summer and fall data, GM_94b outperforms gNDVI significantly with a lower RMSE and considerably higher R 2 (Table 9). This highlights the need for multiple validation dates to properly select the inversion method if one aims to monitor canopies over time.
The number of best matches q to keep when doing the LUT-based inversion also varied depending on the validation date, and a generalizable method to determine it was not found. There is still no clear methodology to determine this number. Banskota et al. [16] used three different values for the inversion (q = 50 when using 6 bands from AVIRIS-C images, q = 150 when using their full spectrum, and q = 200 when using Landsat images; LUT size: 115,000); Atzberger and Richter [64] only considered the best match when using 8 bands from CHRIS images (LUT size: 1,056,000); Weiss et al. [65] used q = 50 when trying to inverse digital data sets (LUT size: 280,000); and Ferreira et al. [66] considered q = 20 when using AISA images (LUT size: 3500). This is unfortunately not satisfactory, as much room is left for subjectivity and overfitting of the inversion method.

Assessment of LAI and Pigment Estimations Accuracy
Even though the site LAI is very low (mean LAI is 0.8 m 2 /m 2 and mean field points' LAI is 0.88 m 2 /m 2 ), LAI estimations were encouraging. The LAI heterogeneity within the Tonzi site was captured, as for a given date, points from both low and high LAI parts of the site were correctly estimated. Seasonal variations were also visible, as LAI variations from summer to fall are clearly followed. Indeed, in Figure 9, blue and dark green points come from the same part of the site (see site map, Figure 1) but are from fall and summer, respectively. Estimation RMSE (0.22 m 2 /m 2 ) is remarkably low, and nRMSE (RMSE divided by the mean field points' LAI) is 0.25, quite similar to what can be obtained over a higher-LAI site: when estimating LAI over a Mediterranean grassland with a mean LAI of 2.87 m 2 /m 2 , Darvishzadeh et al. [67] obtained nRMSEs of 0.22 and 0.18 using statistical and RTM-inversion approaches, respectively.
C ab estimation performances were in line with previous studies. Foliar chlorophylls content estimations (RMSE = 5.94 µg/cm 2 , R 2 = 0.75) are as good as what could be obtained over higher-LAI coniferous and deciduous sites in various studies (RMSE = 4.34 µg/cm 2 , R 2 = 0.47 [68]; RMSE = 4.2 µg/cm 2 [69]; RMSE = 8.1 µg/cm 2 , R 2 = 0.40 [29]). Carotenoid estimations did not perform that well, even though the estimation RMSE was low (RMSE = 2.57 µg/cm 2 , R 2 = 0.1). However, Figure 10b shows that the low R 2 is mostly due to the dark orange point which is, as for C ab , severely underestimated. Further, the foliar Car estimation of the other points appears to be acceptable. Using high-resolution imagery (50 cm), Zarco-Tejada et al. [57] obtained an RMSE below 1.3 µg/cm 2 and R 2 of at most 0.46 when using the SAILH and the FLIGHT radiative transfer models for carotenoid estimation over vineyards. One must also consider that the Car variation range of the present study goes from 5 to 13 µg/cm 2 , while the LUT step is only 4 µg/cm 2 : despite this, the R 2 values obtained are in line with those obtained by Zarco-Tejada et al. [57].
Another factor that could explain the estimation errors (and specifically the underestimation of the dark orange point's biochemistry) is the assumption that the biochemistry measured on one tree can be considered representative of (i) the biochemistry of leaves that are sensed by AVIRIS-C (i.e., mostly top of canopy leaves) and (ii) the leaf biochemical content of all the trees within a pixel. While to address (i), attention was given to only collecting leaves on the sunlit part of the crown (see Section 2.2.2), it is not possible to confirm (ii) with the available field data, and a leap of faith must be made. There is therefore always an uncertainty concerning how representative validation data are with regards to their associated pixels.

Challenges Concerning LAI, C ab and Car Monitoring
Concerning seasonal monitoring, while the QUDO phenophase is such that the site LAI and foliar C ab should have presented a maximum in summer [70,71], Figure 12 showed a decrease from spring to fall. The most probable cause is that the C3 grass was not completely dead in May, as such a behavior is consistent with the results found by Stöckli et al. [72]. This induces a combined contribution of both the green grass and the QUDO vegetation layers in the spring. Even though grass reflectance was taken into account within DART models so that its biochemistry would not influence results, the pure soil pixels manually selected from the AVIRIS-C image for ground reflectance input within the DART model may have contained less green grass than other, non-pure-ground pixels. This would lead to erroneous estimations where green grass is present. As for both summer and fall, the grass is completely dry, which was not an issue for these seasons, and a LAI and foliar C ab decrease was detected. Foliar carotenoid content was correctly found to be increasing between summer and fall for the two years with available fall images, which allows for a secondary validation of carotenoid estimations for this period. Concerning the high foliar Car during spring 2015, a possible cause would be the still-dying grass during the AVIRIS-C acquisition, due to problems similar to those already mentioned for LAI. Estimated LAI and C ab variations were greater than the estimation RMSE found in Sections 3.3.2 and 3.3.3, which gives confidence in the ability of the chosen methodology to monitor the site over summer and fall.
While this study used clear-sky images, clouds are often present in hyperspectral images. Obviously too thick clouds hide the ground, and cloud shadows can be hard to detect and correct for, but models to correct for thin cirrus presence in the 0.4-2.5 µm spectral range exist [73][74][75] and could be used in case of operational use.

Conclusions
Estimating canopy pigments and LAI using radiative transfer inversion is challenging in highly heterogeneous environments. This study investigated the ability of LUT-based inversion of AVIRIS-C medium resolution (18 m) hyperspectral images to estimate LAI and biochemical properties (C ab , Car) of low-foliage open canopies (mean LAI: 0.82 m 2 /m 2 , mean CC: 47%) in a woodland savanna. The methodology followed during this study used almost no a priori knowledge, so as to be able to monitor the study site in its entirety over multiple dates and obtain more general results. Various inversion criteria were tested, relying either on spectral intervals or hyperspectral vegetation indices. Results from very different site locations in terms of LAI, canopy cover, and tree structure were consistent and showed good accuracy for LAI and leaf C ab retrieval and were also encouraging concerning leaf Car retrieval. These results are further strengthened by the temporal aspect of the estimations, with validation done using both summer and fall data from various years. Estimated LAI and C ab also showed a temporal pattern coherent with expected seasonal variations for dates when field data were not available.
This study found that the best inversion criterion for LAI estimation is: NDVI ( ρ 0.833µm −ρ 0.677µm ρ 0.833µm +ρ 0.677µm ) difference; for leaf C ab : GM_94b ( ρ 0.750µm ρ 0.550µm ) difference; and for leaf Car: RMSE between pixel and LUT reflectances over the 0.5-0.-55 µm spectral interval. These criteria showed results either comparable to or below previous studies done over denser canopies (RMSE for LAI: 0.22 m 2 /m 2 ; for C ab : 5.94 µg/cm 2 ; for Car: 2.57 µg/cm 2 ). Surprisingly, the vegetation indices presenting the best performances were not specifically designed to be soil-adjusted, which could indicate that the modeling within DART was sufficiently realistic when it came to ground reflectance, despite the ground being modeled as a lambertian surface with uniform reflectance.
The results demonstrated that medium resolution (18m) hyperspectral imagery, combined with a 3D RTM only using a simple geometric representation of the canopy, is appropriate for the monitoring of vegetation bio-physicochemical parameters for low-LAI open broadleaved forests such as woodland savannas. In the future, the proposed estimation methodology should be further tested on another Californian site (San Joaquin Experimental Range, https://www.neonscience.org/field-sites/fieldsites-map/SJER) to confirm these results. Estimations should also be extended to leaf water content and leaf mass per area so as to have a complete overview of the forest health over time.  Appendix A Figure A1. Estimated canopy LAI from 2013 to 2018. Areas in white were masked as they contained roads, buildings, water, too few QUDO and PISA trees or too many PISA trees. Figure A2. Estimated canopy C ab from 2013 to 2018. Areas in white were masked as they contained roads, buildings, water, too few QUDO and PISA trees or too many PISA trees. Figure A3. Estimated canopy Car from 2013 to 2018. Areas in white were masked as they contained roads, buildings, water, too few QUDO and PISA trees or too many PISA trees.