Sensitivity of Estimated Total Canopy SIF Emission to Remotely Sensed LAI and BRDF Products

International Institute for Earth System Sciences, Nanjing University, Nanjing, Jiangsu, China Yuxiu Postdoctoral Institute, Nanjing University, Nanjing, Jiangsu, China Jiangsu Provincial Key Laboratory of Geographic Information Science and Technology, Key Laboratory for Land Satellite Remote Sensing Applications of Ministry of Natural Resources, School of Geography and Ocean Science, Nanjing University, Nanjing, Jiangsu, China Huangshan Park Ecosystem Observation and Research Station, Ministry of Education, China Department of Geography and Planning, University of Toronto, Toronto, Ontario, Canada Max Planck Institute for Biogeochemistry Department Biogeochemical Integration, Germany

However, only a portion of fluorescence, which is originally emitted by chlorophyll-a molecules in the photosynthesis system [17,18], escapes from canopies and then is observed by sensors (SIF obs ) in a particular observation direction [19][20][21]. The difference in escape probability among biomes could also cause the difference in the GPP-SIF obs relationship. For example, SIF escapes less from needle leaf forest than broadleaf forest canopies due to the higher clumping effect of needle leaf forest [22], and this difference in SIF escape probability between these two types of forests should be considered in the relationship between GPP and SIF obs .
Therefore, it is required to estimate the total canopy SIF emission (SIF total ) to mitigate the canopy structural and angular effects on the estimation of GPP from satellite SIF data [20,[23][24][25][26]. Several studies have proposed methods to derive SIF total from SIF obs using statistically based approaches, such as random forest algorithm [23], and physically based approaches, such as the spectral invariant theory that requires canopy interception (i 0 ) and canopy reflectance of vegetation (R V ) [19,25]. Regardless of statistically or physically based approaches, auxiliary data such as MERIS Terrestrial Chlorophyll Index (MTCI) and canopy reflectance at bands of 685 nm, 710 nm, and 785 nm are required by Liu et al. [23] and near-infrared reflectance of vegetation (NIR V ) and leaf area index (LAI) are required by Zhang et al. [26]. Due to its simplicity and efficiency in deriving SIF total , the approach based on NIR V and i 0 has been adopted by Zhang et al. [26] to derive global SIF total from OCO-2 SIF obs . A more consistent relationship between GPP and SIF total across C 3 plants is established, demonstrating the advantage of SIF total for global GPP estimation.
Although better relationships of SIF total with GPP than SIF obs have been reported for TROPOMI [20] and OCO-2 SIF [26], uncertainties in the above-mentioned satellite products are still considerable, which could impact the relationships between GPP and SIF total . Currently, the trade-off between the advantage of accounting for the escape probability and the disadvantage of the uncertainty in the auxiliary data has not been well investigated to better understand the usefulness of SIF total . Nevertheless, it is difficult to accurately estimate the uncertainties for all satellite products. As a powerful tool, the Soil-Canopy-Observation of Photosynthesis and the Energy balance (SCOPE) model [27] can capture the physical mechanisms behind photosynthesis and fluorescence, and it has been extensively used in the community of SIF remote sensing [28][29][30][31]. Therefore, the SCOPE model can be used to simulate the uncertainty effect on the relationships between GPP and SIF total by artificially adding random uncertainty.
The fluorescence spectrum emitted by chlorophyll-a molecules in the 650-850 nm range has two peaks in the red (~685 nm, RSIF) and far-red (~740 nm, FRSIF) [17,18,32]. Both RSIF and FRSIF originate from photosystem I (PS I) and photosystem II (PS II) [17]. RSIF is mainly from PS II, which is better linked to photochemical quenching and nonphotochemical quenching [33,34]. As expected, RSIF should be more sensitive to GPP than FRSIF [35]. This is supported by a global sensitivity analysis of the SCOPE model [36]. In addition, Zhang et al. [37] also reported that RSIF shows better seasonal correlation with photosynthesis than FRSIF from Scots pine at the leaf scale during the spring recovery of photosynthesis. Furthermore, canopy SIF total at both red (RSIF total ) and far-red band (FRSIF total ) has also been investigated with field observations [24,38], but their performance in estimating GPP has not been inves-tigated and compared with satellite observations. Moreover, the signal of RSIF is weaker than that of FRSIF due to the stronger reabsorption of pigments, which reduces the retrieval accuracy of RSIF compared to FRSIF [9,39].
In this work, the main objectives are (1) to investigate the effect of uncertainties in i 0 and R V on the relationships between GPP and SIF total based on the SCOPE model simulations and (2) to evaluate the sensitivity of estimated SIF total to uncertainties in multiple satellite LAI and BRDF products.

Materials and Methods
2.1. TROPOMI SIF Data. Both RSIF and FRSIF from TRO-POMI were used in this study (ftp://fluo.gps.caltech.edu/ data/). The Sentinel 5 Precursor (S-5P) satellite with a single payload of TROPOMI was launched on 13 October 2017 on a near-polar, sun-synchronous orbit. The repeat cycle in the nadir direction is 17 days, and the overpass time at equator is~13 : 30 local time. S-5P has a varying across track spatial resolutions of 3.5-14 km according to pixel position but a fixed along track spatial resolution of 7.2 km (5.6 km after 6 August 2019). Recently, TROPOMI SIF has been successfully retrieved using a data-driven approach based on a singular value decomposition technique in the atmospheric windows of 663-685.3 nm for RSIF [39] and 743-758 nm for FRSIF [9]. Details about the retrieval process can be referred to Köhler et al. [9] and Köhler et al. [39] and hence not shown here for simplicity.

Derivation of Total
Canopy SIF Emission at the Photosystem Level (SIF total ). A full description of the theoretical basis behind the derivation of SIF total can be found in recent studies [19,20,25,40]. Only a brief description is presented here. The escape probability of fluorescence from leaf surface to canopy (f LC ) for dense canopies and black soil can be approximately estimated as follows [19]: where R λ is the bidirectional reflectance factor in the same wavelength (λ) and observation direction as SIF obs and ω λ is leaf albedo (leaf reflectance + transmittance). To reduce soil effects on R λ for sparse canopies, Zeng et al. [25] proposed to replace R λ at near-infrared band with nearinfrared reflectance of vegetation (NIR V ), which is the product of reflectance in near-infrared band and normalized difference vegetation index (NDVI) [41]: where R nir and R red are the reflectance at near-infrared and red bands, respectively. Similarly, red reflectance of vegetation (Red V ) was calculated with R red and NDVI [24]: Journal of Remote Sensing To avoid confusion, R red is the reflectance of whole canopy (vegetation + soil) in the red band and Red V is the reflectance of vegetation in the red band. Since no corresponding reflectance data is currently available for satellite SIF at the same sun-viewing geometry as SIF, the RossThick-LiSparseR (RTLSR) BRDF model can be used to simulate reflectance at red and near-infrared bands that can be further used to calculate Red V and NIR V . Therefore, Red V and NIR V can be used as R λ in Equation (1) for calculating f LC for RSIF and FRSIF, respectively. Parameters to drive the RTLSR BRDF model can be available from existing BRDF products and details are presented in Section 2.3. i 0 is commonly calculated with G-function (G), leaf area index (LAI), clumping index (CI), and solar zenith angle (SZA, θ) as follows [42]: where the empirical derived parameters ϕ 1 and ϕ 2 are dependent on χ L , which is the departure of leaf angles from a random distribution, and χ L is assigned as biome-specific values based on the Common Land Model 4.5 (CLM 4.5) [43]. We derive the global values of χ L based on MODIS plant functional type classification (MCD12Q1), and the spatial maps of χ L can be found in Figure 1. The CI data was from He et al. [22]. Details of LAI products used in this study are presented in Section 2.4. The sensitivities of the calculation of SIF total to different BRDF and LAI products were systematically evaluated to serve as a reference for the calculation of satellite SIF total . To derive SIF total at the photosystem level, the escape probability of fluorescence from photosystem to the leaf surface (f PL ) was introduced [24]. Therefore, the escape probability of SIF (f PC or f esc ) from photosystem level to canopy level in any direction is calculated as follows: Both ω λ and f PL are negatively related to the absorptance of pigments, such as chlorophylls a + b (Cab). In other words, high Cab causes low ω λ and f PL , and vice versa. To simplify Equation (5), we define K λ as the ratio of ω λ to f PL and assume K λ can be roughly estimated with Cab. Based on SCOPE model simulations (see Section 2.6), K λ for RSIF at 683 nm quickly decreased with Cab and started to saturate when Cab > 40 μg/cm 2 ( Figure 2). In comparison, K λ for FRSIF at 740 nm showed less variations within the ranges of 1.2-1.6. Due to the lack of accurate Cab information, we simply set K λ as 0.6 and 1.2 for RSIF and FRSIF, respectively, which are suitable for a wide range of Cab.

Journal of Remote Sensing
Finally, SIF total at the photosynthesis level can be calculated as follows:

Bidirectional Reflectance Distribution Function (BRDF)
Parameter Products. To be consistent with the same sunviewing geometry as TROPOMI SIF, reflectance at red and far-red bands was simulated using the semiempirical models required for Red V and NIR V calculation, such as RossThick-LiSparseR (RTLSR) as follows: where θ, υ, and ϕ are the solar zenith, view zenith, and relative azimuth angles, respectively. The first term (f iso ) on the right-hand side of Equation (7) represents Lambertian reflectance. f vol and f geo are the coefficients for volumescattering (K vol ) and geometric-optical (K geo ) kernels, respectively. These coefficients (f iso , f vol , and f geo ) are available for three BRDF products, including MCD43A1, VNP43IA1, and MCD19A3, used in this study. Several major information (such as spatial and temporal resolutions) about these products is listed in Table 1, and more details (such as retrieval algorithm) can be found in the listed references. These coefficients provided by both MCD43A1 and VNP43IA1 were derived using the top-ofcanopy reflectance with varying sun-target-viewing geometries after atmospheric correction [44]. The coefficients in MCD19A3 were directly derived from top-of-atmosphere L1B reflectance using the MultiAngle Implementation of Atmospheric Correction (MAIAC) algorithm [45,46]. Both MCD43A1 and VNP43IA1 were released at a daily interval, while MCD19A3 was released in an 8-day interval. For all three products, the BRDF parameters with best quality were used in this study according to the QA layer. The Red V (NIR V ) derived from MCD43 BRDF, VNP43 BRDF, and MCD19 BRDF were denoted as MCD43 Red V (NIR V ), VNP43 Red V (NIR V ), and MCD19 Red V (NIR V ), respectively. In addition, the differences in band configurations between TROPOMI and MODIS/VIIRS sensors were ignored due the marginal RMSE < 0:007 and 0.04 for red and NIR bands, respectively ( Figure 3).

Leaf Area Index (LAI) Products.
Three LAI products were used, including MODIS LAI (MCD15A2H), VIIRS LAI (VNP15A2H), and CGLS LAI (GEOV2) (see details in Table 2). MCD15A2H and VNP15A2H retrieval algorithms are based on a 3-D radiative transfer model that can simulate spectral canopy properties for each biome [47,48]. A lookup-table technique was developed as the main method to retrieve LAI. When the main method failed, a back-up solution based on the empirical relationships between LAI and NDVI was used [48]. Note that only LAI retrievals from the main method were used in this study. CGLS LAI (version GEOV2) was derived from PROBA-V using an artificial neural network (ANN) that was trained based on MODIS/TERRA collection 5 and CYCLOPES V3.1 data [49,50]. The LAI values outside the expected ranges were excluded according to the quality flag (QFLAG) provided in the CGLS products. The temporal series of CGLS LAI were smoothed, with a temporal resolution of 10 days [51], and MCD15A2H and VNP15A2H were composited over 8 days [52]. The uncertainty in i 0 (σ i0 ) was calculated using the error propagation model as follows: where σ LAI is the retrieval uncertainty in LAI products. In this study, σ LAI was obtained from the standard deviation provided in MCD15 and VNP15 LAI products and the RMSE (root mean square error) provided in CGLS LAI product. The i 0 derived from MCD15 LAI, VNP15 LAI, and CGLS LAI were denoted as MCD15 i 0 , VNP15 i 0 , and CGLS i 0 , respectively. All absolute uncertainties in LAI and i 0 were divided by their own values to represent the relative uncertainties (in %) following Fang et al. [53].  Table 3). The standard gap-filling approach was applied to half-hour flux (such as net ecosystem CO 2 exchange) and meteorological data (such as air temperature, vapor pressure deficit, and shortwave incoming radiation) [54]. Subsequently, the gap-filled data were used to calculate half-hourly GPP with the night-time partitioning procedures, in which the daytime respiration was estimated from air temperature using the model calibrated with nighttime data [55]. For each flux site, the monthly TROPOMI SIF was determined as the mean value of all cloud-free observations (cloud fraction < 0:2) within a 10 km radius of the site location. These days with cloudy fraction < 0:2 were denoted as clear-sky days. The half-hour GPP data on clear-sky days were averaged to monthly GPP to match with the satellite SIF. Similarly, LAI and BRDF data for each site were also aggregated to a 10 km radius to be consistent with SIF.
2.6. SCOPE Model Simulation. The effects of the uncertainty on i 0 , NIR V , and Red V for the performance of SIF total in GPP estimation were first analyzed using the SCOPE model (v1.73) [27] before analyzing the satellite SIF data and in situ GPP. The SCOPE model can simulate both SIF and GPP, providing a tool to investigate the relationships between GPP and two SIF metrics (SIF obs and SIF total ). We simulated     Journal of Remote Sensing 5000 scenarios with the random combinations of biochemical, structural, and meteorological parameters listed in Table 4. To be consistent with TROPOMI SIF, the simulated RSIF obs and FRSIF obs were extracted at narrow bands centered at 683 nm and 740 nm, respectively. Different levels of random uncertainty ranging from 0 to 40% were added to i 0 , NIR V , and Red V to investigate the sensitivity of SIF total to the uncertainty in remote sensing products. Note that the SCOPE simulations can be considered the instantaneous observations for GPP and SIF [27]. In addition, only the C 3 photosynthesis pathway was considered for simplicity, because similar results can be expected between C 3 and C 4 photosynthesis pathways. The sun and view geometric information was represented as solar zenith angle, view zenith angle, and relative azimuthal angle. For each simulation scenario, random combinations of all parameters in their own ranges (Table 4) were generated.

Results
3.1. Sensitivity of the GPP-SIF total Relationships to Uncertainty in i 0 , Red V , and NIR V . The relationships of instantaneous GPP with RSIF obs and FRSIF obs based on the SCOPE simulations are shown in Figure 4(a) and 4(b), in which hyperbolic models were suitable for capturing the nonlinearity. Without considering the variation in the escape probability, RSIF obs was weakly and nonlinearly related to GPP (R 2 = 0:38, Figure 4(a)), and FRSIF obs was moderately and nonlinearly related to GPP (R 2 = 0:65, Figure 4(b)). These R 2 for RSIF obs vs. GPP and FRSIF obs vs. GPP were set as the benchmark to evaluate the usefulness of SIF total after considering the escape probability effect. When the uncertainties were not added into Red V , NIR V , and i 0 , both RSIF total and FRSIF total exhibited improved relationships with GPP. R 2 increased from 0.38 for GPP vs. RSI-F obs to 0.76 for GPP vs. RSIF total (Figure 4(c)), and R 2 increased from 0.65 for GPP vs. FRSIF obs to 0.79 for GPP and FRSIF total (Figure 4(d)). The SCOPE simulation demonstrated the usefulness of SIF total to improve the link to GPP by accounting for the varying escape probability. Figure 5 shows the 2-D distribution of R 2 for hyperbolic models between GPP and SIF total derived from i 0 , Red V , and NIR V with different levels of uncertainties. In general, R 2 decreased with the increased level of uncertainties in i 0 , Red V , and NIR V . The black lines in Figure 5 represent the contour lines with R 2 of 0.38 for RSIF obs and 0.65 for FRSI-F obs . As compared to RSIF obs , RSIF total would be well related to GPP when the uncertainty in i 0 and Red V is less thañ 30% ( Figure 5(a)). Similarly, if the uncertainty in i 0 and NIR V is less than~20%, FRSIF total would also be better related to GPP than FRSIF obs (Figure 5(b)). However, when the uncertainties in Red V , NIR V , and i 0 exceeded the uncertainty threshold (~30% for RSIF and~20% for FRSIF), the estimated SIF total was too noisy and could not improve the relationships with GPP compared to SIF obs .

Comparison of i 0 , Red V , and NIR V among Different
Products. High consistencies were found among MCD15 i 0 , VNP15 i 0 , and CGLS i 0 ( Figure 6). The correlation between MCD15 i 0 and VNP15 i 0 was as high as 0.99 (Figure 6(a)),   Incoming shortwave radiation (0.4-2.5 μm) 100-1000 7 Journal of Remote Sensing which was expected due to their similar retrieval algorithms for LAI. In contrast, CGLS i 0 exhibited good but weaker relationships with MCD15 i 0 (R 2 = 0:92 in Figure 6(b)) and VNP15 i 0 (R 2 = 0:90 in Figure 6(c)). In addition, three LAI products showed similar levels of uncertainty after normalizing uncertainties with each LAI itself (Figures 7(a)-7(c)). The uncertainties in percentage were estimated as 23.51%, 22.95%, and 22.86% for MCD15 LAI, VNP15 LAI, and CGLS LAI, respectively (Figures 7(a)-7(c)). Therefore, i 0 derived from these LAI products also showed similar levels of uncertainty in the range of 16.27%-17.10% (Figures 7(d)-7(f)). Fortunately, the uncertainty levels of i 0 derived from all three satellite LAI products were less than the thresholds (30% for RSIF and 20% for FRSIF) determined by the SCOPE simulations ( Figure 5).
Moderately to highly strong relationships were found among Red V derived from three BRDF products with R 2 in the range of 0.82 to 0.93, and only a few data points diverged from the regression line (Figures 8(a)-8(c)). This indicated that these BRDF products estimated consistent vegetation  Journal of Remote Sensing reflectance in the red band overall. Compared to the red band, stronger and more consistent relationships were found for NIR V , with R 2 > 0:95 (Figures 8(d)-8(f)). Since these BRDF products did not provide uncertainty information, the uncertainty of Red V and NIR V was not compared in this study.

Tower Flux GPP against TROPOMI RSIF and FRSIF.
The scatter plots of GPP against RSIF obs and FRSIF obs from TROPOMI are shown in Figure 9, in which C 3 and C 4 plants were separated due to their distinct photosynthe-sis pathways. The nonlinear model was used in the instantaneous GPP and SIF based on the SCOPE simulations, but linear models were efficient to capture the relationships between monthly GPP and SIF for both C 3 and C 4 plants. In general, FRSIF obs showed better relationships with GPP than RSIF obs regardless of C 3 and C 4 plants, which was consistent with the SCOPE simulations. In addition, higher R 2 and slope of linear model were observed for C 4 plants than that for C 3 plants. This higher slope for C 4 plants is attributed to the lower photorespiration and higher efficiency of photosynthesis in plants with  9 Journal of Remote Sensing C 4 metabolism than C 3 plants [56,57]. We also observed an interesting phenomenon that the intercept was positive and negative for C 3 and C 4 plants, respectively. Theoretically, both SIF and GPP are from APAR, so the intercept should be zero (when APAR is zero) under natural conditions. The nonzero intercept reported here could be caused by the regression model uncertainties, the bias in satellite SIF retrievals, the bias in flux tower GPP partition, and the environmental stress.
After accounting for the difference in escape probability, the relationships of GPP with RSIF total and FRSIF total are presented in Figures 10 and 11. The poorest relationship was obtained between GPP and RSIF total in C 3 plants with R 2 from 0.55 to 0.57 (Figure 10), which was still higher than the R 2 between GPP and SIF obs (R 2 = 0:49 in Figure 9(a)). Subsequently, R 2 between GPP and FRSIF total in C 3 plants ranged between 0.72 and 0.77 (Figure 11), which outperformed the relationship between GPP and FRSIF obs 10 Journal of Remote Sensing (R 2 = 0:70 in Figure 9(b)). As for C 4 plants, RSIF total also improved the relationships with GPP from RSIF obs ( Figure 10). However, FRSIF total did not show improvement in R 2 for C 4 plants (Figure 11), since the R 2 between FRSI-F obs and GPP has reached to 0.88 (Figure 9(b)). Although there was no clear difference in R 2 obtained by different combinations of i 0 and Red V or NIR V , we observed slightly higher R 2 obtained by MCD19 Red V and NIR V than those by MCD43 and VNP43 in terms of the relationship between GPP and FRSIF total . Figure 10: Relationships between tower-based GPP and TROPOMI RSIF total calculated with different combinations of Red V and i 0 . Red V was calculated with MCD43 (first row), VNP43 (second row), and MCD19 (third row). i 0 was calculated with MCD15 (first column), VNP15 (second column), and CGLS (third column). For simplicity, these sites were classified into two types: C 3 and C 4 plants, with following number indicating R 2 . All regression models are statistically significant (p value < 0.001).

Journal of Remote Sensing
The escape probability of SIF (fesc) estimated with Equation (5) is also compared with that from the SCOPE simulations ( Figure 12). For RSIF at 685 nm, fesc calcu-lated with Red V and i 0 from different combinations of BRDF and LAI products was clearly higher than that from the SCOPE simulations. Therefore, RSIF total (=RSIF obs /fesc) MCD19 MCD19 Figure 11: Relationships between tower-based GPP and TROPOMI FRSIF total calculated with different combinations of NIR V and i 0 . NIR V was calculated with MCD43 (first row), VNP43 (second row), and MCD19 (third row). i 0 was calculated with MCD15 (first column), VNP15 (second column), and CGLS (third column). For simplicity, these sites were classified into two types: C 3 and C 4 plants, with following number indicating R 2 . All regression models are statistically significant (p value < 0.001).

12
Journal of Remote Sensing was underestimated in this study. A further work is still required to reduce the underestimation of RSIF total . In contrast, fesc at 740 nm from SCOPE simulations was better consistent with that calculated with NIR V and i 0 , demonstrating the success of NIR V + i 0 in calculating fesc at 740 nm.

Comparison between the SCOPE Simulation and
TROPOMI SIF total . The superiority of SIF total in GPP estimation has been shown by recent studies [20,23,24]. However, the sensitivity of SIF total to uncertainty in i 0 , NIR V , and Red V has not been well understood. Based on the SCOPE simulations, the improvement in R 2 (ΔR 2 ) can reach up to 0.38 and 0.14 for RSIF total and FRSIF total , respectively, when no uncertainty existing in i 0 , NIR V , and Red V (Figure 4). However, the differences in R 2 tend to decrease with the increasing level of uncertainty in i 0 , NIR V , and Red V as shown in Figure 5, revealing the adverse effect of uncertainty on the relationships between GPP and SIF total . Since the uncertainty in satellite data is unavailable, ΔR 2 for actual scenarios is likely to be less than the maximum values used in this study. For example, the actual ΔR 2 only ranges between 0.02 (2.86%) and 0.07 (10.00%) for TROPOMI FRSIF total in C 3 plants (Figure 11), which is close to the reported ΔR 2 (0.04, 5.40%) obtained by OCO-2 FRSIF total [26]. The lower ΔR 2 in actual scenarios is attributed to the uncertainty in i 0 , NIR V , and Red V . The uncertainties from clumping index, Gfunction, and leaf albedo can also contribute to the lower ΔR 2 , although LAI and BRDF products are the main source of the uncertainties in SIF total . The comparison between simulation and measurement promotes our understanding of the use of SIF total . Furthermore, this study discusses the potential uncertainty in LAI and BRDF products as below.

Impacts of Different
Satellite LAI and BRDF Products on the Estimation of SIF total . Numerous studies have intercompared existing LAI products from regional to global scales in terms of spatiotemporal consistency and reported many difference among these products [58][59][60][61][62][63][64][65][66][67][68]. These inherent uncertainties in LAI products result from both retrieval algorithms and input data [48]. For example, to improve the inversion efficiency for MCD15 and VNP15, several biome-specific variables (e.g., canopy structure, leaf type, and soil brightness) are defined beforehand in the inversion process. As a result, the biome-specific assignments could result in uncertainty for LAI retrievals for mixed or misclassified pixels [48]. In addition, the uncertainty in atmosphere parameters could propagate into the atmospheric correction process [69], bringing also uncertainty to reflectance and hence LAI retrievals. For most studies, CGLS product shows better accuracy as compared to MCD15 and VNP15. For example, Brown et al. [58] reported the better agreements between reference LAI and CGLS LAI than MCD15 and VNP15 LAI. However, this study observes consistent relationships among MCD15, VNP15, and CGLS i 0 with R 2 ranging from 0.90 to 0.99 ( Figure 6) and similar uncertainty level (~17%) in i 0 ( Figure 7). As expected, improved estimation of GPP from SIF total is available if uncertainty in i 0 is further reduced. This can be obtained by new satellite sensors with improved spectral and spatial resolutions and more accurate retrieval algorithm, such as ESA's forthcoming FLuorescence EXplorer (FLEX) mission in tandem with Sentinel-3 [70]. Since i 0 and fraction of absorbed photosynthetically active radiation (FAPAR) are highly related [25], SIF total calculated with FPAR also exhibits similar results as compared SIF total calculated with i 0 (results not shown).
Several studies also reported the high consistency between MCD43 and VNP43 NDVI [71] and MCD43 and MCD19 NIR [72], which supports the high consistency in NIR V (the product of NDVI and NIR reflectance) among MCD43, VNP43, and MCD19 ( Figure 8). Therefore, marginal differences are expected in relationship between GPP and SIF total calculated from different Red V and NIR V (Figures 10 and 11). In terms of FRSIF total , the slightly higher R 2 for MCD19 than those for MCD43 and VNP43 could be attributed to the advantage of MAIAC algorithm adopted by MCD19. In addition, MODIS is onboard two satellites (Terra and Aqua) with two equator local crossing times of 10:30 and 13:30, and VIIRS is only onboard one satellite (Suomi-NPP) with an equator local crossing time of 13:30; the former could provide more angular samplings and full inversion for BRDF parameters than the latter [73]. However, the equator local crossing time for VIIRS is consistent with that for TROPOMI SIF, which should be more suitable for TROPOMI SIF than MCD43 and MCD19. With the additional VIIRS launched in 2017 and to be launched in the future as part of the JPSS program, an increased pixel number of full inversions will be available to generate the VNP43 product. As a result, reducing uncertainty in VNP43 could improve the calculation of SIF total and the estimation of GPP.

Conclusions
Previous studies have shown that SIF total was more useful for GPP estimation than SIF obs across multiple scales. However, the advantage of SIF total in improving GPP estimation could be masked by the uncertainty in the derivation of i 0 , NIR V , and Red V , which were required by the calculation of SIF total . In this study, we first investigated the effect of the uncertainty in i 0 , Red V , and NIR V on the calculation of SIF total and the relationships between SIF total and GPP based on the SCOPE model simulations. As a result, SIF total performed better than SIF obs for both red and far-red bands in capturing the link with GPP. The improvement in R 2 (ΔR 2 ) for SIF total and GPP relationships was 0.38 and 0.14 for RSIFtotal and FRSIF total from RSIF obs and FRSIF obs , respectively. With the increasing uncertainty in i 0 , NIR V , and Red V , RSIFtotal and FRSIF total showed degraded relationships with GPP. Furthermore, ΔR 2 decreased to zero when the uncertainty levels were higher than~30% in i 0 and Red V (for estimation of RSIF total ) and~20% in i 0 and NIR V (for estimation of FRSIF total ) based on the SCOPE model simulation. Then, this study calculated RSIF total and FRSIF total from TROPOMI RSIF obs and FRSIF obs with different combinations of i 0 (from MCD15, VNP15, and CGLS LAI) and Red V and NIR V (from MCD43, MCD19, and VNP43). In general, TROPOMI RSIFtotal and FRSIF total exhibited better relationships with flux tower GPP than RSIF obs and FRSIF obs . Due to the comparable uncertainty levels among these different satellite products (such as LAI), the estimation of SIF total was less sensitive to the choice of satellite products. Our results based on SCOPE simulations and TROPOMI data contribute to our understanding of the estimation of SIF total using current satellite products, which would advance the use of satellite SIF data for global terrestrial GPP estimation.