Retrieval and validation of forest background reflectivity from daily Moderate Resolution Imaging Spectroradiometer (MODIS) bidirectional reflectance distribution function (BRDF) data across European forests

Biogeosciences, 18, 621–635, 2021 https://doi.org/10.5194/bg-18-621-2021 © Author(s) 2021. This work is distributed under the Creative Commons Attribution 4.0 License. Retrieval and validation of forest background reflectivity from daily Moderate Resolution Imaging Spectroradiometer (MODIS) bidirectional reflectance distribution function (BRDF) data across European forests Jan Pisek1, Angela Erb2, Lauri Korhonen3, Tobias Biermann4, Arnaud Carrara5, Edoardo Cremonese6, Matthias Cuntz7, Silvano Fares8, Giacomo Gerosa9, Thomas Grünwald10, Niklas Hase11, Michal Heliasz4, Andreas Ibrom12, Alexander Knohl13, Johannes Kobler14, Bart Kruijt15, Holger Lange16, Leena Leppänen17, Jean-Marc Limousin18, Francisco Ramon Lopez Serrano19, Denis Loustau20, Petr Lukeš21, Lars Lundin22, Riccardo Marzuoli9, Meelis Mölder4, Leonardo Montagnani23,31, Johan Neirynck24, Matthias Peichl25, Corinna Rebmann11, Eva Rubio19, Margarida Santos-Reis26, Crystal Schaaf2, Marius Schmidt27, Guillaume Simioni28, Kamel Soudani29, and Caroline Vincke30 1Tartu Observatory, University of Tartu, Tõravere, Tartumaa, Estonia 2School for the Environment, University of Massachusetts Boston, Boston, Massachusetts, USA 3School of Forest Sciences, University of Eastern Finland, Joensuu, Finland 4Lund University, Lund, Sweden 5Fundación CEAM, Paterna, Valencia, Spain 6ARPA Valle d’Aosta, Saint-Christophe, Italy 7Université de Lorraine, AgroParisTech, INRAE, UMR Silva, Nancy, France 8CNR – National Research Council, Rome, Italy 9Department of Mathematics and Physics, Università Cattolica del Sacro Cuore, Brescia, Italy 10Institute of Hydrology and Meteorology, Department of Hydro Sciences, Technische Universität Dresden, Dresden, Germany 11Helmholtz Centre for Environmental Research – UFZ, Leipzig, Germany 12Department of Environmental Engineering, Technical University of Denmark, Kongens Lyngby, Denmark 13Faculty of Forest Sciences and Forest Ecology, University of Göttingen, Göttingen, Germany 14Umweltbundesamt GmbH, Vienna, Austria 15Department of Environmental Sciences, Wageningen University & Research, Wageningen, the Netherlands 16Norwegian Institute of Bioeconomy Research, Ås, Norway 17Space and Earth Observation Centre, Finnish Meteorological Institute, Sodankylä, Finland 18CEFE, Université Montpellier, CNRS, EPHE, IRD, Université Paul-Valéry Montpellier, Montpellier, France 19IER-ETSIAM, Universidad de Castilla-La Mancha, Albacete, Spain 20INRAE, Bordeaux, France 21Global Change Research Institute, Academy of Sciences of the Czech Republic, Brno, Czech Republic 22Department of Soil and Environment, Swedish University of Agricultural Sciences, Uppsala, Sweden 23Faculty of Science and Technology, Free University of Bolzano, Bolzano, Italy 24INBO, Geraardsbergen, Belgium 25Department of Forest Ecology and Management, Swedish University of Agricultural Sciences, Umeå, Sweden 26cE3c – Centre for Ecology, Evolution and Environmental Changes, Lisbon, Portugal 27Forschungszentrum Jülich, Jülich, Germany 28INRAE URFM, Avignon, France 29Université Paris-Saclay, CNRS, AgroParisTech, Ecologie, Systématique et Evolution, Orsay, France 30Faculty of Bioscience Engineering, Earth and Life Institute, Université Catholique de Louvain, Louvain-la-Neuve, Belgium Published by Copernicus Publications on behalf of the European Geosciences Union. 622 J. Pisek et al.: Retrieval and validation of forest background reflectivity 31Forest Services, Autonomous Province of Bolzano, Bolzano, Italy Correspondence: Jan Pisek ([email protected]) Received: 30 September 2020 – Discussion started: 4 November 2020 Revised: 12 January 2021 – Accepted: 12 January 2021 – Published: 27 January 2021 Abstract. Information about forest background reflectance is needed for accurate biophysical parameter retrieval from forest canopies (overstory) with remote sensing. Separating underand overstory signals would enable more accurate modeling of forest carbon and energy fluxes. We retrieved values of the normalized difference vegetation index (NDVI) of the forest understory with the multi-angular Moderate Resolution Imaging Spectroradiometer (MODIS) bidirectional reflectance distribution function (BRDF)/albedo data (gridded 500m daily Collection 6 product), using a method originally developed for boreal forests. The forest floor background reflectance estimates from the MODIS data were compared with in situ understory reflectance measurements carried out at an extensive set of forest ecosystem experimental sites across Europe. The reflectance estimates from MODIS data were, hence, tested across diverse forest conditions and phenological phases during the growing season to examine their applicability for ecosystems other than boreal forests. Here we report that the method can deliver good retrievals, especially over different forest types with open canopies (low foliage cover). The performance of the method was found to be limited over forests with closed canopies (high foliage cover), where the signal from understory becomes too attenuated. The spatial heterogeneity of individual field sites and the limitations and documented quality of the MODIS BRDF product are shown to be important for the correct assessment and validation of the retrievals obtained with remote sensing. 1 Introduction The reflectance from the forest canopy background/forest floor can often confound and even dominate the radiometric signal from the upper forest canopy layer to the atmosphere. Forest understory is defined here as all the components found under the forest canopy, including understory vegetation, leaf litter, moss, lichen, rock, soil, snow, or a mixture thereof (Pisek and Chen, 2009). If unaccounted for, forest understory can introduce potential bias in the estimation of overstory biophysical parameters (e.g., leaf area index, LAI, and fraction of absorbed photosynthetically active radiation, fAPAR) and, subsequently, productivity estimates (e.g., the net primary productivity – NPP) as the contribution of the understory to the total energy absorption capacity of a forest stand can be quite significant (Clark et al., 2001; Law et al., 2001). The understory vegetation in forest ecosystems should be treated differently from the overstory in carbon cycle modeling because of the different residence times of carbon fixed through NPP in different ecosystem components (Vogel and Gower, 1998; Rentch et al., 2003; Marques and Oliveira, 2004; Kim et al., 2016). Currently, the understory is often treated as an unknown quantity in carbon models due to the difficulties in measuring it properly and consistently across larger scales (Luyssaert et al., 2007). The predictions regarding the spectral variation in forest background have posed a persistent challenge (McDonald et al., 1998; Gemmell, 2000) because of the high variability in incoming radiance below the forest canopy, challenges with the spectral characterization, and weak signal in some parts of the spectrum for both overstory and understory (Schaepman et al., 2009), and the general varying nature of the understory (Miller et al., 1997). Multi-angle remote sensing can capture signals of different forest layers because the observed proportions for different forest layers vary with the viewing angle, making it possible to separate forest overstory and understory signal. Here, we aim to consolidate previous efforts of tracking understory reflectance and its dynamics with multi-angle Earth observation data (Canisius and Chen, 2007; Pisek and Chen, 2009; Pisek et al., 2010, 2012, 2015a, 2015b, 2016; Jiao et al., 2014) by testing the validity of this approach, using Moderate Resolution Imaging Spectroradiometer bidirectional reflectance distribution function (MODIS BRDF)/albedo data (gridded 500m daily Collection; 6 MCD43 product), against in situ understory reflectance measurements over an extended set of Integrated Carbon Observation System (ICOS) forest ecosystem sites. The validation procedure was defined to comply as much as possible with the best practices proposed by the Committee on Earth Observation Satellites (CEOS) Working Group on Calibration and Validation (WGCV) Land Product Validation (LPV) subgroup (Garrigues et al., 2008; Baret et al., 2006). It corresponds to Stage 1 validation, as defined by the CEOS (Nightingale et al., 2011; Weiss et al., 2014), where product accuracy shall be assessed over a small (typically <30) set of locations and time periods by comparison with in situ or other suitable reference data. Using the extended set of ICOS forest ecosystems as validation sites, we asked the following questions: 1. Can ICOS forest ecosystem sites serve as a suitable validation data set with respect to their footprint and the pixel resolution of Earth observation (EO) products? Biogeosciences, 18, 621–635, 2021 https://doi.org/10.5194/bg-18-621-2021 J. Pisek et al.: Retrieval and validation of forest background reflectivity 623 Figure 1. Distribution of study sites across Europe; for further details, refer to Table 1. 2. Can we retrieve reliable normalized difference vegetation index (NDVI; Rouse et al., 1973; Tucker, 1979) dynamics for understory with MODIS BRDF data across diverse forest conditions during the growing season? 3. Are there important differences between the total (overstory and understory) and understory-only NDVI signals? 2 Materials and methods 2.1 Study sites The ICOS is a distributed pan-European research infrastructure providing in situ standardized, integrated, long-term and high-precision observations of lower atmosphere greenhouse gas (GHG) concentrations and land–atmosphere and ocean–atmosphere GHG interactions (Gielen et al., 2017). In this study, we carried out the evaluation over the network of 31 ICOS-affiliated forest ecosystem sites, complemented with additional sites in Spain, Portugal, Austria, and Finland. Together, these selected 40 study sites comprise a large variety of forest overand understory types, spanning a wide latitudinal gradient from almost 38◦N (Yeste, Spain) to 68◦N (Kenttärova, Finland). Site locations are shown in Fig. 1, and vegetation characteristics are summarized in Table 1. 2.2 Understory spectra and forest canopy cover/closure in situ measurements Following the terminology by Schaepman-Strub et al. (2006), we refer to the reflectance factors measured by the field spectrometers as the satellite-derived hemispherical–directional reflectance factors (HDRFs). The given spectrometer’s field of view is approximated as being angular (cone) and narrower than a whole hemisphere, with some anisotropy captured which corresponds to normal remote sensing viewing geometry. An overview of the undertaken in situ campaigns at each site and their characteristics are given in Table 1. The individual sites were visited between April 2016 and August 2019, mostly during the growing season. Following the protocol by Rautiainen et al. (2011), the understory spectra were measured with the Sun completely obscured by the clouds or at around sunset (diffuse light conditions), covering the visible/near infrared (NIR) region, depending on the spectrometer (see Table 1 for more details). A total of three understory spectra were measured every 2m along two 50m long transects laid at each site, resulting in 50 measurement points (150 individual measurements). Transects covered and characterized conditions within the measurement footprint of the given tower. It should be noted that the tower footprint might be different from the exact MODIS pixel footprint (see Sect. 2.5 for the spatial homogeneity assessment of MODIS pixels). The measurements concerned conditions on the forest floor and low herbaceous and shrubby species or tree seedlings and saplings, as the area sampled by each spectral measurement was estimated to correspond to a ∼50cm diameter circle on the ground. The downward-pointing spectroradiometer (no fore optics were used) was held by the operator’s outstretched hand. A total of three spectra above a 10in. (0.254m) Spectralon SRT-99–100 white panel were recorded at the beginning, after every four understory spectra measurement points (every 8m), and at end of each transect. A hemispherical, conical reflectance factor was obtained with an uncalibrated Spectralon reflectance spectrum and the linearly interpolated irradiance. Finally, broadband HDRFs for red (620–670nm) and NIR (841–876nm) wavelengths were computed with relative spectral response functions for the MODIS sensor onboard Terra. The understory NDVI value for given site was calculated from the red and NIR-band values and averaged over the two transects. Estimates of overstory foliage cover and crown cover were obtained from digital cover photographs (DCPs). Overstory foliage cover was defined as the percentage of ground covered by the vertical projection of foliage and branches and crown cover as the percentage of ground covered by the vertical projections of the outermost perimeters of the crowns on the horizontal plane (without double-counting the overlap; Gschwantner et al., 2009). The DCPs were taken from below the canopy every 8m along transects at each site. The camera (Nikon CoolPix4500; 2272×1704 resolution) was https://doi.org/10.5194/bg-18-621-2021 Biogeosciences, 18, 621–635, 2021 624 J. Pisek et al.: Retrieval and validation of forest background reflectivity Table 1. Study site characteristics and their spatial representativeness status. ICOS – Integrated Carbon Observing System sites; LTER – Long Term Ecological Research Network sites. Note that the sampling dates are shown in the format yyyy/mm/dd. Site code Site name Lat (◦) Long (◦) Sampling date Spectrometer model Understory vegetation Representativeness AT-Zbn Zöbelboden (LTER) 47.842 14.442 2017/11/18 ASD FieldSpec 4 Calamagrostis varia, Brachypodium sylvaticum, Hordelymus europaeus, and Senecio ovatus Not representative BE-Bra Brasschaat (ICOS) 51.304 4.519 2019/01/12 Ocean Optics; FLAME-S-VIS-NIR-ES Betula spec., Quercus robur, and Sorbus aucuparia Representative BE-Vie Vielsalm (ICOS) 50.3 5.983 2018/08/16 Ocean Optics; FLAME-S-VIS-NIR-ES Sparse fern and moss cover Representative at 0.5km CH-Dav Davos (ICOS) 46.817 9.85 2018/07/12 Ocean Optics; FLAME-S-VIS-NIR-ES Dwarf shrubs, blueberry, and mosses Representative CZ-BK1 Bílý Kˇ ríž (ICOS) 49.502 18.539 2016/04/17 ASD FieldSpec 4 Vaccinium myrtillus L. Representative at 1.5km CZ-Lnz Lanžhot (ICOS) 48.682 16.948 2017/04/27 ASD FieldSpec 4 Allium ursinum and Asarum europeum Representative DE-Hai Hainich (ICOS) 51.079 10.453 2018/04/12 Ocean Optics; FLAME-S-VIS-NIR-ES Anemone nemorosa and Allium ursinum Representative DE-HoH Hohes Holz (ICOS) 52.083 11.217 2018/04/11 Ocean Optics; FLAME-S-VIS-NIR-ES Anemone nemorosa Representative DE-RuW Wüstebach (ICOS) 50.505 6.331 2018/08/16 Ocean Optics; FLAME-S-VIS-NIR-ES Sparse Deschampsia flexuosa, Deschampsia cespitosa and Molinia caerulea Not representative DE-Tha Tharandt (ICOS) 50.967 13.567 2018/04/12 Ocean Optics; FLAME-S-VIS-NIR-ES Fagus sylvatica,Abies alba, and Deschampsia flexuosa Representative DK-Sor Soroe (ICOS) 55.486 11.645 2018/09/26 Ocean Optics; FLAME-S-VIS-NIR-ES Beech saplings and seedlings; Pteridium aquilinum Representative ES-AP1 Almodóvar del Pinar 39.677 −1.848 2017/11/09 ASD FieldSpec HandHeld 2 Quercus ilex ssp. ballota,Rosmarinus officinalis, Thymus vulgaris, Lavandula latifolia, Quercus coccifera, and Genista scorpius Representative ES-CMu Cuenca del Majadas 40.252 −1.965 2017/11/12 ASD FieldSpec HandHeld 2 Juniperus communis,Juniperus oxycedrus, and Crataegus monogyna Representative at 0.5km ES-CPa Cortes de Pallas 39.224 −0.903 2017/11/08 Ocean Optics; FLAME-S-VIS-NIR-ES Rosmarinus officinalis,Ulex parviflorus, and Brachypodium retusum Representative >0.5km ES-Yst Yeste 38.339 −2.351 2018/07/28 Ocean Optics; FLAME-S-VIS-NIR-ES Rosmarinus officinalis L.,Thymus vulgaris L., and Cistus clusii Dunal Representative at 0.5km FI-Hal Halssiaapa 67.368 26.654 2017/06/13 ASD FieldSpec Pro Sedge vegetation Representative FI-Hyy Hyytiälä (ICOS) 61.847 24.295 2018/06/28 Ocean Optics; FLAME-S-VIS-NIR-ES Vaccinium spec. and Norway spruce seedlings Representative at 0.5km Biogeosciences, 18, 621–635, 2021 https://doi.org/10.5194/bg-18-621-2021 J. Pisek et al.: Retrieval and validation of forest background reflectivity 625 Table 1. Continued. Site Code Site Name Lat (◦) Long (◦) Sampling date Spectrometer model Understory vegetation Representativeness FI-Ken Kenttärova (ICOS) 67.987 24.243 2017/06/13 ASD FieldSpec Pro Vaccinium myrtillus,Empetrum nigrum, and Vaccinium vitis-idaea and the forest mosses Pleurozium schreberi, Hylocomium splendens, and Dicranum polysetum Representative at 2km FI-Kns Kalevansuo 60.647 24.356 2017/06/15 ASD FieldSpec Pro Dwarf shrubs and mosses Representative at 0.275km FI-Let Lettosuo (ICOS) 60.642 23.96 2017/06/15 ASD FieldSpec Pro Dwarf shrubs, mosses, and herbs Representative <1.0km FI-Sod Sodankylä (ICOS) 67.362 26.638 2017/06/13 ASD FieldSpec Pro Lingonberry, Calluna vulgaris, and lichens Spheroid does not fit <1.5km; not representative at >1.5km FI-Var Värriö (ICOS) 67.757 29.616 2017/06/14 ASD FieldSpec Pro Mosses, lichens, and dwarf shrubs Representative FR-Bil Bilos – Salles (ICOS) 44.494 −0.956 2018/06/14 Ocean Optics; FLAME-S-VIS-NIR-ES Molinia coerulea Moench.,Pteridium aquilineum, and Ulex europaeus Representative <0.5km FR-FBn Font Blanche (ICOS) 43.241 5.679 2018/06/12 Ocean Optics; FLAME-S-VIS-NIR-ES Quercus coccifera,Phillyrea latifolia, and other species Representative at 1.5km FR-Fon Fontainebleau–Barbeau (ICOS) 48.476 2.780 2018/06/16 Ocean Optics; FLAME-S-VIS-NIR-ES Carpinus betulus Representative at 0.5km FR-Hes Hesse (ICOS) 48.674 7.066 2018/08/18 Ocean Optics; FLAME-S-VIS-NIR-ES Fagus sylvatica seedlings and blackberry Spheroid does not fit FR-MsS Montiers (ICOS) 48.537 5.312 2019/01/14 Ocean Optics; FLAME-S-VIS-NIR-ES Sparse Sphagnum spec. vegetation Representative >1.0km FR-Pue Puéchabon (ICOS) 43.741 3.596 2018/06/13 Ocean Optics; FLAME-S-VIS-NIR-ES Buxus sempervirens,Pistacia lentiscus, Phillyrea latifolia,Salvia rosmarinus, and Ruscus aculeatus Representative IT-BFt Bosco Fontana (ICOS) 45.202 10.743 2018/07/10 Ocean Optics; FLAME-S-VIS-NIR-ES Hedera helix,Corylus spec., and Ruscus aculeatus Representative at 1.5km IT-Cp2 Castelporziano 2 (ICOS) 41.704 12.357 2019/01/25 Ocean Optics; FLAME-S-VIS-NIR-ES Phyllirea latifolia and Pistacia lentiscus Representative at 1.5km IT-Ren Renon (ICOS) 43.732 10.291 2018/07/11 Ocean Optics; FLAME-S-VIS-NIR-ES Deschampsia flexuosa L., Vaccinium myrtillus L., and Rhododendron ferrugineum L. Representative at 0.5km IT-SR2 San Rossore (ICOS) 61.847 24.295 2018/06/28 Ocean Optics; FLAME-S-VIS-NIR-ES Ligustrum vulgare Representative <1.5km IT-Trf Torgnon 45.833 7.567 2018/07/07 Ocean Optics; FLAME-S-VIS-NIR-ES Juniperus communis,Rhododendron ferrugineum, and Festuca varia Not representative https://doi.org/10.5194/bg-18-621-2021 Biogeosciences, 18, 621–635, 2021 626 J. Pisek et al.: Retrieval and validation of forest background reflectivity Table 1. Continued. Site Code Site Name Lat (◦) Long (◦) Sampling date Spectrometer model Understory vegetation Representativeness NL-Loo Loobos (ICOS) 52.167 5.744 2018/08/13 Ocean Optics; FLAME-S-VIS-NIR-ES Prunus serotina,Vaccinium myrtillus, Deschampsia flexuosa, and mosses Representative at 0.5km NO-Hur Hurdal (ICOS) 60.372 11.078 2018/09/27 Ocean Optics; FLAME-S-VIS-NIR-ES Vaccinium spec. and Norway spruce seedlings Representative PT-Cor Coruche (LTER) 39.138 -8.333 2016/10/08 Ocean Optics; FLAME-S-VIS-NIR-ES Rumex acetosella,Tuberaria guttata, Tolpis barbata Plantago coronopus, Agrostis pourretii,Briza maxima, Vulpia bromoides, and Vulpia geniculata Spheroid does not fit SE-Htm Hyltemossa (ICOS) 56.098 13.419 2018/09/28 Ocean Optics; FLAME-S-VIS-NIR-ES Continuous moss cover Spheroid does not fit SE-Knd Kindla (LTER) 59.754 14.908 2016/07/16 ASD FieldSpec Pro Ericaceous dwarf shrubs, mosses, and lichens Representative SE-Nor Norunda (ICOS) 60.086 17.48 2018/10/22 ASD FieldSpec Pro Bilberry, lingonberry, and moss Representative <1.5km SE-Svb Svartberget (ICOS) 64.256 19.775 2019/08/23 ASD FieldSpec Pro Bilberry, lingonberry, and moss Representative set to automatic exposure, aperture-priority mode, minimum aperture, and F2 lens (Macfarlane et al., 2007). The camera was leveled at the height of 1.4m above the ground, and the lens was pointed towards the zenith. This setup provides a view zenith angle from 0 to 15◦, which is comparable with the first ring of the LAI-2000 instrument (Macfarlane et al., 2007). We used the algorithm by Nobis and Hunziker (2005) to threshold the majority of the DCP images. However, some of the images were visibly overexposed, i.e., the 8bit digital numbers (DNs) of the background sky were 255, and parts of any portion of the sky were black (typically at 240–250 DN). Next, a method based on mathematical image morphology (Korhonen and Heikkinen, 2009) was applied to estimate the foliage and crown cover fractions. In this method, black and white canopy images are processed with morphological closing and opening operations that are well known in digital image processing (Gonzalez and Woods, 2002). As a result, a filter for large gaps was obtained. When a tuning parameter (called the structuring element in image processing) was set so that large gaps only occurred between individual tree crowns (Korhonen and Heikkinen, 2009), the proportions of gaps inside and between individual crowns could be calculated. 2.3 Background signal retrieval method with EO data The total reflectance of a pixel (R) results from the weighted linear combination of reflectance values by the forest canopy, forest background, and their sunlit and shaded components (Li and Strahler, 1985; Chen et al., 2000; Bacour and Bréon, 2005; Chopping et al., 2008; Roujean et al., 1992) as follows: R=kTRT+kGRG+kZT RZT +kZGRZG,(1) which includes the reflectivities of the sunlit crowns (RT), sunlit understory (RG), shaded crowns (RZT ), and shaded understory (RZG). RGmarks the bidirectional reflectance factor (BRF) of the target (understory). The kjare the proportions of these components at the chosen viewing angle or in the instantaneous field of view of the sensor at the given irradiation geometry. Following Canisius and Chen (2007), we derive the understory reflectivity (RG) with the assumption that the reflectivities of the overstory and understory at the given illumination geometry differ little between the chosen viewing angles. While the components may not fully meet the definition of Lambertian reflectors (i.e., reflecting electromagnetic radiation equally in all directions), several previous studies (e.g., Bacour and Bréon, 2005; Deering et al., 1999; Peltoniemi et al., 2005) found forward-scattering reflectance factors of various targets off the principal plane to be fairly constant. The most suitable viewing configuration for the retrieval has been identified by Pisek et al. (2015a), using a high angular resolution BRF data set of Kuusk et al. (2014) and accompanying in situ measurements of understory reflectance factors (Kuusk et al., 2013). The configuBiogeosciences, 18, 621–635, 2021 https://doi.org/10.5194/bg-18-621-2021 J. Pisek et al.: Retrieval and validation of forest background reflectivity 627 ration consists of the BRF at nadir (Rn=0◦), with the solar zenith angle (SZA) corresponding to the Sun’s position at 10:00 local time (LT) for given day and another zenith angle (Ra=40◦) with relative azimuth angle PHI=130◦. It can be expressed by Eqs. (2) and (3) as follows: Rn=kT nRT+kGnRG+kZT nRZT +kZGnRZG (2) Ra=kT aRT+kGaRG+kZT aRZT +kZGaRZG.(3) The proportions of the components (kj) were obtained using the four-scale model (Chen and Leblanc, 1997) with parameters for generalized deciduous and coniferous tree stands as an input (see Table 2; Kuusk et al., 2013). The understory reflectance at the desired wavelengths can be calculated by combining and solving Eqs. (2) and (3) and the insertion of Rnand Raestimates derived from appropriate EO data. The individual components (sunlit/shaded overstory and understory) cannot be resolved with the MODIS spatial resolution. The reflectances of shaded tree crowns (RZT ) and understory (RZG) are related to sunlit ones via Mas RZT =M×RTand RZG =M×RG, where M=RZ/R for a reference target, which can be measured in the field or predetermined with the four-scale model. Here, the same Mis assumed for overstory trees and the understory. Based on his field work in Canadian boreal forests, White (1999) suggested that angularly constant, wavelength-dependent Mvalues may be appropriate, at least during the growing season. The input stand parameters from Table 2 may not be always precisely known while retrieving the understory signal over larger areas. Figure 2 shows the relationships between the available in situ data for tree heights or tree densities over our study sites with the 1km2resolution estimates from the global maps of Simard et al. (2011) and Crowther et al. (2015). The weak relationships indicate the current unsuitability of the site-specific variable estimates of interest (tree height and tree density) from currently available global maps at a given spatial resolution for our purpose. At the same time, the calculated mean values for the tree heights of needleleaved (17.5m) and broadleaved tree stands (22.7m) from Simard et al. (2011) over the study sites were reasonably close to our original generalized input parameter values in Table 2. Following Gemmell (2000), we opted to report a range of understory NDVI (NDVIu) values obtained with the combination of parameter values from Table 2 for each site and date. Specifying the correct constraints (window) for background alone has been previously found to greatly reduce the errors in the estimation of overstory parameters (Gemmell, 2000). 2.4 MODIS BRDF data The MCD43A1 V6 bidirectional reflectance distribution function and albedo (BRDF/albedo) model parameter data set is a 500m gridded daily product. MCD43A1 is generated by inverting multi-date, multi-angular, cloud-free, atmospherically corrected, and surface reflectance observations acquired by MODIS instruments onboard the Terra and Aqua Figure 2. (a) Relationship between available in situ estimates of tree height (in meters) with Simard et al.’s (2011) estimate. (b) Relationship between available in situ estimates of tree density (trees per hectare) with Crowther et al.’s (2015) estimates. DBF – deciduous broadleaf forest; EBF – evergreen broadleaf forest; DNF – deciduous needleleaf forest; ENF – evergreen needleleaf forest; MF – mixed forest. satellites over a 16d period (Wang et al., 2018). The Julian date represents the ninth day of the 16d retrieval period, and consequently, the observations are further weighted to estimate the BRDF/albedo for that particular day of interest. The MCD43A1 algorithm uses all high-quality observations that adequately sample the viewing hemisphere to fit an appropriate semiempirical BRDF model (the RossThickLiSparseReciprocal model; Roujean et al., 1992; Lucht et al., 2000) for that location and date of interest. We computed the bidirectional reflectance factor (BRF) at the top of the canopy with the isotropic parameter and two (volumetric and geometric) kernel functions (Roujean et al., 1992) for MODIS band 1 (red – 620–670nm) and band 2 (NIR – 841–876 nm). We used the Ross and Li kernels to reconstruct the BRF values for required geometries (see Sect. 2.3) for each date, and then we derived the understory signal, using the formulas described in Sect. 2.3. The associated data quality (MCD43A2) product was employed to assess the effect of the retrieval quality on the accuracy of the calculated understory signal. All MODIS data have been accessed and processed through the Google Earth Engine (Gorelick et al., 2017). 2.5 Spatial representativeness assessment of the validation sites A method developed by Román et al. (2009), and refined by Wang et al. (2012, 2014, 2017) was adopted to evaluate the spatial representativeness of in situ measurements to assess the uncertainties arising from a direct comparison between field-measured forest understory spectra and the corresponding estimates with MODIS BRDF data. To characterize the spatial representativeness of a test site to represent a satellite retrieval, this method uses three variogram model parameters (the range, sill, and nugget), obtained by the analysis of near-nadir surface reflectances from cloud-free 30m Landsat/Operational Land Imager (OLI) data (Román et al., 2009) collected as close to the sampling date as possible. https://doi.org/10.5194/bg-18-621-2021 Biogeosciences, 18, 621–635, 2021 628 J. Pisek et al.: Retrieval and validation of forest background reflectivity Table 2. Stand parameters for the four-scale model. Stand Deciduous Coniferous Stand density (treesha−1) 500, 1000, and 2000 500, 1000, and 2000 Tree height (m) 25 16 Length of live crown (m) 9.2 4.2 Radius of crown projection (m) 1.87 1.5 Leaf area index (m2m−2) 1, 2, and 3 1, 2, and 3 Figure 3. Shortwave BRF composites centered at ICOS sites of (a) Norunda in Sweden and (c) Wüstebach in Germany. (b, d) Variogram estimators (points), spherical model results (dotted curves), and sample variances (solid straight lines) obtained over the sites with Operational Land Imager (OLI) subsets and spatial elements of 0.275, 0.5, 1.0, 1.5, and 2km as a function of the distance between observations. Variogram legend explanations: a – variogram range; var – sample variance; c – variogram sill; c0 – nugget variance. Where valid imagery was not available within a reasonable window of the sampling date, imagery from the corresponding season of a different year was used. As such, the analysis was done to illustrate the representativeness of the tower site with respect to a particular point in time. Campagnolo et al. (2016) showed that the effective spatial resolution of 500m gridded MODIS BRDF product at mid-latitudes is around 833×618m because of the varied footprints of the source multi-angular surface reflectance observations. We analyzed each site with five different spatial extents (0.275, 0.5, 1, 1.5, and 2km) to assess and illustrate the changes in spatial representativeness with different spatial resolutions. 3 Results and discussion 3.1 Spatial representativeness Table 1 provides the assessment of spatial heterogeneity for all sites included in this study, using OLI subsets acquired around the time of in situ measurements. The example results for the ICOS sites of Norunda (SE-Nor) in Sweden and Biogeosciences, 18, 621–635, 2021 https://doi.org/10.5194/bg-18-621-2021 J. Pisek et al.: Retrieval and validation of forest background reflectivity 629 Wüstebach (DE-RuW) in Germany, using three OLI subsets, are shown in Fig. 3. The variogram functions with relevant model parameters for the two sites are displayed in Fig. 3b and d. The range corresponds to the value on the xaxis where the model flattens out. There is no further correlation of a biophysical property associated with that point beyond the range value. The sill is the ordinate value of the range. A smaller sill value indicates a more homogenous surface (less variation in surface reflectance). A surface can be considered spatially representative with respect to the MODIS footprint when the sill value is <5.0e−4(Román et al., 2009; Wang et al., 2017). The sill values for all spatial extents are well below the value of 5.0e−4, up to 1km spatial resolution in the case of Norunda (Fig. 3b), which indicates that the field measurements are representative and allow comparison with MODIS retrievals at a 500m spatial resolution. While the Wüstebach site can be considered spatially homogeneous within the immediate vicinity of 275m around the tower, the sill value exceeds the criteria of 5.0e−4 at >0.5km spatial resolution. During late summer/early autumn of 2013, trees were almost completely removed in an area of 9ha west of the tower in order to promote the natural regeneration of a near-natural deciduous forest from a spruce monoculture forest. The clearfelling area can be seen in Fig. 3d. This action resulted in an increase in the spatial heterogeneity of this ICOS site. In situ measurements collected within the footprint of the Wüstebach tower, thus, cannot be deemed fully comparable with the retrievals with MODIS at a 500m spatial resolution. Overall, most of the sites were found representative at the spatial resolution of MODIS BRDF gridded data. The nonrepresentative cases and the effect on the understory signal retrieval and agreement with the corresponding in situ measurements carried within the measurement footprint of the individual towers are further discussed in Sect. 3.2 and 3.3. Román et al. (2009) provide further details on the assessment of spatial representativeness, using a set of four geostatistical attributes derived from semivariograms. 3.2 NDVI ranges There is only a weak relationship between the total (overstory and understory) NDVI signal retrieved with MODIS BRDF data and corresponding in situ understory NDVI measurements (R2=0.19; Fig. 4). Total NDVI values alone do not allow one to disentangle the correct understory signal. In contrast, our retrieval method could track the understory signal dynamics over a broad NDVI range (Fig. 5). The predicted understory NDVI ranges were beyond the uncertainty limits of in situ understory measurements (corresponding to ±1 standard deviation (SD) here) in less than 15% of cases. These sites with poor retrievals were carefully investigated to identify the issues precluding good results. Below, we focus on a discussion of results where the predicted and in situ measured NDVI ranges of the understory layer did not agree. Figure 4. Relationship between total (overstory and understory) NDVI values computed from nadir NDVI values, using MODIS BRDF/albedo data and in situ measured understory NDVI values over the study sites. The understory dominated the overall signal of open shrubland at the Cortes de Pallas (ES-CPa) site and the deciduous broadleaf forest site at Montiers (FR-MsS) during the leaf-off part of the season (Fig. 5). Both sites were found to be spatially representative for comparison with MODIS footprint data at the time of the available in situ measurements (Table 1). There are only very few trees scattered across the Cortes de Pallas site, and ground vegetation is fully exposed. Extremely low tree density does not match with any of the original generalized input parameter values in Table 2, and the predicted understory signal does not match well with the in situ measurements. In situ measurements at Montiers were carried out during the leaf-off part of the season, which allowed a full exposure of the understory. Despite this, the predicted understory NDVI range from the MODIS data did not overlap with the in situ measurements at Montiers at all. However, the MODIS BRDF values for these sites were marked with lower data quality flags (QA >1), which correctly signals a decrease in accuracy in the calculations of the understory reflectance as well. Overall, our results confirm that, under conditions of very low tree density/leaf-off conditions, the understory signal can be assumed to be identical to the total scene NDVI. The performance of the method turns out to be limited over sites with a closed canopy, such as Bílý Kˇ ríž (CZBK1), Hesse (FR-Hes), or Vielsalm (BE-Vie; Fig. 5). This is because the shadowing effect makes diffuse scattering the dominant mechanism in such stands, and the understory carries only a negligible influence on the top-of-canopy signal. Bosco Fontana (IT-BFt) is another broadleaf forest site with very high foliage cover (FC=0.91), yet the predicted understory NDVI range entirely overlaps with the collected in situ values. It should be noted that, in contrast to other sites with https://doi.org/10.5194/bg-18-621-2021 Biogeosciences, 18, 621–635, 2021