The Effect of Diffuse Radiation on Ecosystem Carbon Fluxes Across China From FLUXNET Forest Observations

Aerosol loading and cloud cover can alter the composition of radiation reaching the Earth’s surface and affect the ecosystem’s carbon cycle. In this study, we established an empirical model of the diffuse radiation fraction (Kd) based on a clearness index (Kt) to obtain the Kd of four FLUXNET forest sites in China. We focused on the relationships among the Kd, photosynthetically active radiation (PAR), light-use efficiency (LUE) and gross primary productivity (GPP) through mechanistic analysis. The relationships between carbon fluxes [including GPP, ecosystem respiration (ER), and net ecosystem exchange (NEE)] and the Kd were explored. Furthermore, we investigated the influence of environmental factors on carbon fluxes. The results showed that the Kd models were accurate in estimating Kd (R2= 0.88–0.93). Overall, the GPP first increased and then decreased with increasing Kd. When Kd< Ko (Ko, the diffuse radiation fraction corresponding to the maximum value of GPP), the direct PAR decreased as Kd increased, while the diffuse PAR increased rapidly. At this stage, the diffuse fertilization effect led to an increase in GPP. When Ko<Kd<Kdiff-max (Kdiff-max, the diffuse radiation fraction corresponding to the maximum value of diffuse PAR), as Kd increased the direct PAR still decreased and the diffuse PAR still increased, but the GPP declined. When Kd>Kdiff-max, the diffuse PAR began to decrease, and the reduction in the superimposed direct PAR caused the GPP of the canopy to drop rapidly. The LUE of the vegetation canopy was higher under diffuse light conditions than under direct light. Furthermore, with an increase in the Kd, the negative value change of the NEE was consistent with the GPP, but the ER was less affected by the Kd. Finally, the impact of temperature (TA) and vapor pressure deficit (VPD) on the GPP was unimodal, and the impact on the NEE was U-shaped. In addition, latent heat (LE) had a significant positive effect on GPP and NEE. Our study emphasized the relationship between the change in PAR composition and the Kd, as well as its impact on the carbon fluxes change, which is highly important to the study of carbon neutralization.


INTRODUCTION
Solar radiation is the main source of energy exchange on the Earth's surface and an important factor affecting the productivity of ecosystems. Among them, wavelengths between 400-700 nm are considered photosynthetically active radiation (PAR), which is the main energy source of plant photosynthesis (Tang et al., 2013;Ren et al., 2018). According to whether they can receive solar radiation, the leaves of the vegetation canopy can be divided into sunlit leaves and shaded leaves . For a canopy with a complex vertical structure, sunlit leaves receive both direct and diffuse radiation and photosynthesis is light-saturated, while shaded leaves receive only diffuse radiation and photosynthesis is usually limited by light (Yue and Unger, 2017;Zhou et al., 2021b;Xue et al., 2021). Aerosol pollution and the presence of clouds increase diffuse radiation (Xue et al., 2020;Zeng et al., 2021). Compared with direct radiation, diffuse radiation has a stronger ability to penetrate the vegetation canopy and a more homogeneous distribution (Ryu et al., 2019), as well as possesses a higher ratio of blue light (Urban et al., 2012), which can promote photosynthesis of shaded leaves and improve the light-use efficiency (LUE), thus improving the gross primary productivity (GPP) of the ecosystem, which is called the diffuse fertilization effect (DFE) (Wang et al., 2015;Wang et al., 2018;Gao et al., 2018;Rap et al., 2018;Xie et al., 2020;Gui et al., 2021). At the same time, however, aerosols and clouds reduce total radiation, and the balance between the drop in GPP caused by the reduction of total radiation and the effect of diffuse fertilization is uncertain (Xing et al., 2017;Han et al., 2019;Zhou et al., 2021a;Xie et al., 2021), which is bad for environment production and economic development in China Yang et al., 2022).
Terrestrial ecosystems absorb carbon dioxide (CO 2 ) through photosynthesis and release carbon into the atmosphere through respiration and therefore play an important role in the global carbon cycle (Xiao et al., 2019). There are many parameters to describe carbon cycle, among which GPP, ecosystem respiration (ER) and net ecosystem exchange (NEE) are the primary component of the land-atmospheric carbon cycle and sequestration (Gui et al., 2021). These biological processes affecting the carbon cycle are not only affected by PAR and diffuse fertilization effects, but also by environmental factors (Normand et al., 2021;Shen et al., 2022a). The efficiency of photosynthetically-transpiration (respiration) is influenced by many factors, including air temperature, solar radiation, water availability, soil temperature, and forest canopy characteristics Shen et al., 2022b). Oliphant et al. (2011) reported that the net effect of vapor pressure deficit (VPD) and air temperature (TA) related to sky conditions contributes to enhancing the light response, although their effect is lower than that of total PAR. Whereas, Dass et al. (2016) showed that the summer temperature explained 37.7% of the annual changes in GPP, while precipitation and cloud accounted for 20.7% and 19.3%, respectively in northern Eurasia. Numerous studies have shown that the effect of environmental factors on GPP cannot be ignored. The change of ecosystem carbon cycle is not the result of the single influence of diffuse radiation, but the result of the combined action of multiple environmental factors, and this influence has great uncertainty, so further exploration is required.
The diffuse radiation fraction (K d ), which is the ratio of diffuse radiation to total radiation, is an important index to measure cloud cover and aerosol load in the sky. It is a key variable to study the influence of diffuse radiation on carbon fluxes and a key process to analyze the impact of PAR and LUE on GPP (Yang et al., 2019;Zhang Y. et al., 2017). Due to limited observations, empirical parameters were historically used to simulate K d (Rensheng et al., 2004;Zeng et al., 2020). The earliest K d model was established by Liu and Jordan based on the clearness index (K t ), which is defined as the ratio of total solar radiation (H) to extraterrestrial radiation (H 0 ) (Liu and Jordan., 1960). Since then, many researchers have used different meteorological parameters and other variables for modeling (Li et al., 2011;Dervishi and Mahdavi, 2012;Cao et al., 2017;Gao et al., 2018;Zhou et al., 2019). According to their different input parameters, these models can be divided into 1) models based on K t , 2) models based on S, 3) models based on K t and S coupling, and 4) models based on multiple parameters (Li et al., 2011;Li et al., 2011). The previous models mostly used quadratic fitting, for instance, Feng et al. (2018) established a variety of functions to simulate K d in mainland China, among which quadratic functions performed best, with the average RMSE, MAE and R were 1.81 MJ m −2 ·d −1 , 1.33 MJ m −2 ·d −1 and 0.79, respectively. But then Yang et al. (2020) compared 18 K d models established in China, and the results showed that the estimated performance of the cubic function model based on K t is better than quadratic function, RMSE is 0.08 MJ m −2 ·d −1 lower. The best-performing model established by Yang et al. (2020) is a binary cubic model coupled with two parameters K t and S. Moreover, the estimated performance of K d model coupled with temperature and relative humidity on the basis of K t and S is slightly better than that of only based on K t , with RMSE 0.06 MJ m −2 ·d −1 lower. But compared with K d model based on K t and S (uncoupled environmental factors), RMSE is 0.02 MJ m −2 ·d −1 higher. In other words, models coupled with environmental factors on the basis of K t and S do not always perform best. Since FLUXNET sites do not provide the observed sunshine duration data and relative humidity data, the K d model of the coupling of S and environmental factors is not suitable for our study, therefore, this study adopted the K d model established on the basis of Kt. With the development of technology, some scholars also use machine learning methods to derive K d (Zhou et al., 2021c). Compared with machine learning methods, the empirical model can take regional radiation conditions and sky conditions into account in the model, which can explain the mechanism of K d to a certain extent, that is, the influence of radiation components and atmospheric transparency on K d . In conclusion, the cubic function model based on K t that we adopted is most appropriate, and can capture the variation of K d better. In addition, compared with machine learning methods, it can explain the mechanism affecting the formation of K d more clearly, and is simpler and more efficient. Moreover, few studies have applied this model to simulate the K d of FLUXNET sites in the past. We tried to build K d model at FLUXNET sites and got good prediction results, which provides a case of model application where there is no diffuse radiation observation data.
The current research on the impact of diffuse radiation on ecosystem productivity has mostly focused on a large spatial area; few studies are based on the mechanisms of observation sites, and few empirical models are built based on FLUXNET sites. Moreover, the relationship among the PAR, GPP and K d is less clear. Whether the increase in the K d is sufficient to compensate for the loss in total PAR and whether the canopy GPP is higher or lower under diffuse radiation need to be further explored. To strengthen research in this field, our study is based on the scale of the site and aims to explore the mechanisms of the interactions between the K d and the PAR, LUE and GPP depending on the sky conditions, vegetation types and environmental conditions of different sites. In this study, the main objectives are to 1) derive the K d value of FLUXNET sites; 2) clarify the complex relationship among K d , PAR, and GPP, as well as K d , LUE, and GPP, and explore the mechanism of the diffuse fertilization effect and the process of GPP change; 3) explore the response of ecosystem carbon fluxes to K d and the influence of other environmental factors on GPP. Forest ecosystems are the main terrestrial carbon sinks, making up approximately 80% of the terrestrial carbon sink (Fang et al., 2018;Piao et al., 2022). The study of the effects of diffuse radiation on forest ecosystem carbon fluxes will provide new insight on carbon peaks and carbon neutrality.

FLUXNET Data
FLUXNET is a global network composed of local regional networks that measures the exchanges of carbon dioxide, water vapor, and energy between the biosphere and atmosphere based on the eddy covariance method (Baldocchi et al., 2001;Friend et al., 2007). The daily GPP, NEE, ER, and heat fluxes, as well as meteorological data (including air temperature, shortwave radiation and VPD) were download from the FLUXNET database (https://fluxnet.org/). There are 12 FLUXNET sites in China, the distribution and overview of which can be seen in Figure 1 and Table 1, respectively. Due to the complex vertical structures of the forest canopy, it has a strong response to diffuse radiation (Strada et al., 2015;Zhang Y. et al., 2020). Therefore, we chose only four forest sites in this study: Cha, Qia, Din and MPM. We extracted flux data from the growing season (June-August) when physiological and phenological changes can be largely neglected (Alton et al., 2007). A detailed description of the geographical environment of the four sites is shown in Table 2. There are three kinds of plant functional types (PFT) at FLUXNET forest sites in China, which are mixed forests (MF), evergreen broadleaf forests (EBF) and evergreen needleleaf forests (ENF). The PFT of different sites is related to the local climate. The Cha site located in a higher latitude zone and has a temperate continental climate, forming a mixed forest of temperate coniferous forest and broad-leaved forest. The other three sites are located in the middle and low latitudes zone and have a subtropical monsoon climate, forming evergreen broad-leaved forests and evergreen coniferous forests. The vertical structure of the canopy plays an important role in diffuse fertilization effect. The response of carbon absorption to diffuse radiation is also related to canopy characteristics such as leaf area index (LAI), leaf angle and leaf optical characteristics. Some studies showed that diffuse fertilization effect increased with the increase of LAI (Gui et al., 2021). From the Table 2, Cha site has the largest LAI and canopy height (5.8 m 2 m −2 and 26m, respectively), followed by Din and Qia sites. MPM sites had the smallest LAI and canopy height (1.75-2.455.8 m 2 m −2 and 2-4 m, respectively). In terms of LAI and canopy height, Cha site has the most complex vertical canopy structure, while MPM site had the simplest canopy structure. In addition, it is worth paying attention to that the typical vegetation type of MPM site is different from other sites. The MPM site is a coastal mangrove ecological community, and mangrove is a woody wetland plant, which has higher requirements on hydrothermal conditions. In addition, the four sites provide data in different years. The Cha, Qia, and Din sites provide data from 2003 to 2005, and the MPM site provides data from 2016 to 2018.

Radiation Data
The in situ-measured radiation dataset was obtained from the National Meteorological Information Center of the China Meteorological Administration (CMA) (https://data.cma.cn/). It provides daily solar radiation data including daily total radiation, diffuse radiation, and sunshine duration. Quality control was carried out for all data to reduce the impact of outliers on the result analysis. Among all these sites, only 19 sites provided diffuse radiation measurements (He and Wang, 2020). The spatial distribution of radiation sites is shown in Figure 1.
For this study, we use daily measurements from 19 diffuse radiation sites during 2000-2017 to establish and verify an empirical model of K d . Among them, the measurements from 2000 to 2011 were used for K d modeling, and data from 2012 to 2017 were used for K d validation. On the other hand, some of these radiation sites are adjacent to FLUXNET sites but they rarely overlap and they are a certain distance away from FLUXNET sites. Therefore, when we built the K d model of FLUXNET sites, we chose the three nearest radiation sites.

Development of K d Model
According to the previous empirical model, we selected the K d model based on K t . K t represents the diffusion and absorption effects caused by the presence of the atmosphere (Bortolini et al., 2013;Khorasanizadeh and Mohammadi, 2016). As K t decreases, sky conditions vary from clear to cloudy (Tsubo and Walker, 2004;Apeh et al., 2021). The FLUXNET sites did not provide diffuse radiation data, but there were many CMA sites around the FLUXNET sites. According to Tobler's First Law of Geography, everything is related to everything else, but things that are located near each other are more related than things located further away (Tobler, 1970). Therefore, we considered the FLUXNET sites as center points and chose the three CMA radiation sites nearest to the FLUXNET sites. The sky conditions and radiation conditions at the radiation sites were used to approximately represent those of FLUXNET sites. We used the data provided by the radiation sites to build the K d model and then input the data from the FLUXNET sites to obtain the corresponding K d . Therefore, according to the distribution of the forest FLUXNET sites, seven CMA sites were chosen. To obtain K d , we first needed   Zhang et al. (2006). For MPM, site information, refer to Wong and Fung (2013).
Frontiers in Earth Science | www.frontiersin.org June 2022 | Volume 10 | Article 906408 to calculate K t , which is the ratio of total solar radiation (H) to extraterrestrial radiation (H 0 ). The equation of daily extraterrestrial radiation (H 0 ) was given as follows (Duffie and Beckman, 2020;Yang et al., 2020): where G SC is the solar constant (1367 W m −2 ), ∅ is the latitude (°), δ is the solar declination (°), ω s is the sunset hour angle (°), and n is the day number of the year. Based on this, we established a regression model of K d and K t in the form of a cubic polynomial function.
Because the measurement of radiation may be subject to errors and failures, it is necessary to perform quality control on the data to identify incorrect data (Yang et al., 2020). This study adopts the data quality control method proposed by Younes (Younes et al., 2005): First, data were rejected if the solar altitude was less than 7°; Second, data were rejected if the ratio of the sunshine duration to maximum possible sunshine was greater than 1; Third, the values of the K t and K d ranged between 0-1; Fourth, the envelope was drawn. The results are shown in Figure 2. After the quality control process, the percentage of valid data at the CMA sites was 95.36-97.34%. The remaining datasets were divided into two parts. The data from 2000 to 2011 were used for K d modeling, and data from 2012 to 2017 were used for K d validation.

Calculation of PAR and LUE
PAR is the main energy source for plant photosynthesis, controlling the rate of effective photosynthesis of terrestrial organisms and directly affecting the growth of vegetation (Li et al., 1997). Its diffuse component can enhance the LUE of the vegetation canopy, thereby increasing carbon uptake. However, not all FLUXNET sites provided PAR observations; only the MPM site provided photosynthesis photon flux density (PPFD), so we converted quantum systems (μmol m −2 s −1 ) into energy systems (W·m −2 ) based on Formula (4). Other sites adopted the method proposed by (Zhu et al., 2010), which was based on geographical division.
where U PAR is the PAR of quantum systems, Q PAR is the PAR of energy systems, and μ is the conversion coefficient, usually 4.55 (Zhou et al., 1984). PAR consists of direct PAR and diffuse PAR, among which diffuse PAR is the most important component that promotes photosynthesis. It is generally estimated by multiplying PAR by the K d (Spitters et al., 1986;Ren et al., 2014). LUE connects solar radiation, photosynthesis and vegetation production, represents the ratio of plant fixation to light energy absorption, and represents the efficiency of converting the light energy captured by plants into chemical energy, that is, the efficiency of CO 2 fixation into organic matter, which can indicate vegetation productivity and carbon fixation efficiency Frontiers in Earth Science | www.frontiersin.org June 2022 | Volume 10 | Article 906408 (Monteith, 1972;Monteith, 1977;Baldocchi and Penuelas, 2019). We attempted to quantify the LUE to analyze the impact of diffuse radiation on the vegetation LUE, which can be defined as the ratio of the GPP to PAR, as shown in Formula (5) (Strada et al., 2015;Yue and Unger, 2017;Wang et al., 2018) (Strada et al., 2015;Yue and Unger, 2017;Wang et al., 2018).

Verification of K d and PAR
To explore the applicability of the K d model at different FLUXNET sites, we evaluated the K d model using observational data from the CMA radiation sites. Due to the proximity of the Din and MPM sites, the three closest CMA radiation sites were the same, so the data used for modeling and validation were the same. The K d model based on the four FLUXNET forest sites is shown in Table 3. Figure 3 shows the density scatterplots of the daily K d simulations from 2012 to 2017. The results illustrated that the K d model had excellent performance in simulating K d , with high R 2 values of 0.88-0.93, low RMSE values of 0.07-0.1 W m −2 , and low MAE values of 0.06-0.07 W m −2 at the four FLUXNET sites. However, at the Cha site, the model overestimated the low value and underestimated the high value of the K d , while at the Qia, Din and MPM sites, there was a slight underestimation of the K d . The land cover, cloud fraction and aerosol conditions varied from site to site, so the model behaved differently across different regions. Compared with the quadratic function model based on K d established by Yang et al. (2020), our RMSE decreased by 0.46-0.81 MJ m −2 ·d −1 , the estimated performance is obviously better.
The PAR data from all the sites except for the MPM site were based on previous empirical formulas. Therefore, it was necessary to simulate PAR. We evaluated the empirical model of PAR using observations at the MPM site because only this site had observation data. Figure 4 shows the density scatterplots of the simulated PAR. The R 2 of the results obtained by the empirical model and the measured results at the FLUXNET site was 0.79, indicating that the model can estimate the PAR of the FLUXNET site well. Its RMSE and MAE were relatively low, at 13.60 and 10.79 W m −2 , respectively. Even if the factors affecting PAR are complex and changeable, the empirical model was based on certain climatological methods, so it could reasonably simulate PAR.

Spatiotemporal Changes in K d
K d is affected by solar radiation, aerosols and clouds, so it changes over time. Exploring the temporal change of the K d helps assess its impact on GPP. We analyzed the differences in FLUXNET sites from the monthly changes in K d . Figure 5 is the monthly variation in the K d simulated at the four forest sites, and the changes in the K d at the four sites are different. At a monthly scale, the K d of the Cha sites was the highest in July and the lowest in January, showing a single-peak trend throughout the year. This is because the Cha site has a typical temperate continental climate, with strong radiation and high cloud cover in summer and the opposite in winter. In contrast to the Cha sites, the Qia sites had the lowest K d in July due to the low cloud cover during the dry season and the highest K d in January, which may be due to the presence of more anthropogenic aerosols in winter. The Din sites had the lowest K d in October due to the low cloud cover during the dry season and the highest K d in  March due to radiation increasing gradually in spring. However, the K d of the MPM site was the highest in April and the lowest in December, but there was little change throughout the year, which was mainly due to the low latitude of the MPM site and the slight annual change in solar radiation. In this study, we mainly studied the effects of diffuse radiation on vegetation productivity during the growing season (summer, including June, July and August), and the K d of the four forest sites in summer from high to low were Din (0.75), MPM (0.66), Qia (0.65) and Cha (0.65). Figure 6 shows the spatial distribution of the K d simulated at all FLUXNET sites. Spatially, the K d of the four seasons was roughly zonal in distribution: high at low latitudes, low at high latitudes, and gradually decreasing from low latitudes to high latitudes. The Tibetan Plateau was also a low value area. This was mainly due to the zonal distribution of solar radiation and cloud cover. Solar radiation in low latitudes is strong, and cloud cover is high, so the diffuse radiation is greater and the K d is higher. The opposite pattern is true at higher latitudes. The Tibetan Plateau has a high altitude, low cloud coverage and strong direct radiation, so the K d is also low. Our research focused on forest sites, but to understand the spatial distribution of the K d at all FLUXNET sites and all vegetation types, it is helpful to understand the K d levels at the four forest sites. Figure 6 shows that the K d values of these four forest sites are higher than those of the other sites in all seasons.

The Relationship Among K d , PAR, and GPP
Changes in aerosol load and cloud cover not only change the PAR reaching the surface but also change the proportion of diffuse radiation, with uncertain overall effects on global plant productivity and the terrestrial carbon sink.
We first explored the relationship among the K d , PAR and GPP, and the result is shown in Figure 7. Clouds and aerosols can increase diffuse radiation, thereby promoting K d . As the K d increased, GPP first increased and then decreased. It was closely related to PAR, and changes in the K d caused changes in PAR and its components. Overall, the total PAR decreased with increasing K d . However, the changes in its direct and diffuse components were more complicated during the different ranges of K d . When K d gradually increased to K o (K d corresponds to the maximum value of GPP), the direct PAR was higher than the diffuse PAR but decreased rapidly, while the diffuse PAR showed  Frontiers in Earth Science | www.frontiersin.org June 2022 | Volume 10 | Article 906408 the opposite results and increased rapidly, and the two radiation components reached a balanced state near K o . At this stage, the growth in GPP was mainly caused by the increase in diffuse PAR. The presence of clouds and aerosols increased the diffuse PAR, and the diffuse PAR in the canopy had a stronger transmittance, which can spread more incident radiation to shaded leaves, greatly promoting the photosynthesis of the shaded leaves. This diffuse fertilization effect leads to an increase in canopy GPP. Due to differences in vegetation PFTs, aerosol load and cloud cover, different regions had an optimal K d (K o ) when the GPP reached a maximum value. The optimal K d (K o ) values of the Cha (MF), Qia (ENF), Din (EBF), and MPM (EBF) sites were 0.54, 0.56, 0.57, and 0.43, respectively, and the sky conditions at this time were the most favorable for carbon assimilation. This result is consistent with the optimal K d range of forests obtained by previous studies, which was between 0.4-0.6 (Alton et al., 2007;Alton, 2008;Zhang J. et al., 2020). When the K d exceeds K o , the amount of aerosols and clouds in the sky is too high, which greatly reduces the total PAR reaching the surface, so the GPP begins to decrease with an increase in K d . However, the root cause of the reduction in total PAR is the substantial reduction in direct PAR. The diffuse PAR first increases and then decreases with increasing K d . This is because the moderate increase in aerosols and cloud cover promotes diffuse PAR. Once the increase becomes excessive, it greatly reduces the total PAR, resulting in no excess PAR becoming scattered. When the K d is between K o and K diff-max (K d corresponds to the maximum value of diffuse PAR), although the diffuse PAR is higher than the direct PAR, it grows slowly, while the direct PAR still decreases, and the diffuse fertilization effect is not enough to compensate for the loss of GPP caused by the reduction of direct PAR, so the overall GPP at this stage has a downward trend. When the K d is greater than the K diff-max , the diffuse PAR also begins to decrease, and because the total PAR has been drastically reduced, there is not enough radiation for the aerosols and clouds to scatter. The drop in diffuse PAR causes a weakening of photosynthesis in shaded leaves, while the direct PAR is still in a reduced state, which leads to a weakening of photosynthesis in sunlit leaves. The reduction in the superposition of the two radiation components results in the rapid decline of canopy GPP.
Hence, the net effect on the canopy carbon assimilation depends on a balance between the reduction in direct PAR (which tends to reduce photosynthesis) and the increase in diffuse PAR (which tends to increase photosynthesis). Figure 7 shows that the LUE increases as the K d increases. Furthermore, we compared the LUE under direct and diffuse radiation conditions. The sky conditions were separated into categories of predominantly diffuse (K d ≥ 0.5) and predominantly direct (K d < 0.5) sunlight. We divided the sky condition based on the K d , and considered direct radiation to be dominant when the K d was less than 0.5 and diffuse radiation to be dominant when the K d was greater than 0.5 (Alton, 2008). The LUE of FLUXNET forest sites under different radiation conditions is shown in Figure 8. The LUE of the four forest sites from high to low was: Cha, MPM, Qia and Din. At all sites, the LUE was significantly higher under diffuse radiation dominated conditions than direct radiation dominated conditions. Specifically, when the sky radiance was dominated by diffuse sunlight, the LUE was found to be 5.52%, 3.38%, 2.58%, and 3.77% higher than that under direct sunlight at the Cha, Qia, Din and MPM sites, respectively. This observation indicated an enhanced LUE, with the greatest increase occurring in dense foliage, such as the Cha site dominated by MF. Table 2 shows that the Cha site has the most complex vertical structure among the four forest sites because it has the highest LAI (5.8 m 2 ·m −2 ) and the highest canopy height (26 m), so it is more sensitive to the diffuse radiation response. Under clear-sky conditions, there were fewer clouds and aerosols, and the K d was low, so solar radiation was dominated by direct light and sunlight was too concentrated on sunlit leaves. This situation led to the sunlit leaves in the canopy being usually light-saturated, and since extra sunlight cannot be used efficiently for photosynthesis, the LUE was low. At the same time, shaded leaves are more efficient at using light, but they were limited by the low incident radiation that they have received. In contrast, under cloudy or aerosol sky conditions, sunlight is more scattered and incoming radiation is more diffuse, the irradiance of the canopy is more uniform, and shaded leaves can receive much more sunlight from all directions, so the lightsaturated part of the canopy is reduced and the LUE increases. This also demonstrates the diffuse fertilization effect. In conclusion, canopy photosynthesis tends to be significantly more light-use efficient under diffuse sunlight than under direct sunlight. The LUE enhancement under diffuse sunlight can be explained by the sharing of the canopy radiation load, which is reduced under direct sunlight. A modest reduction in direct radiation does not reduce photosynthesis in sunlit leaves  Frontiers in Earth Science | www.frontiersin.org June 2022 | Volume 10 | Article 906408 9 because the LUE remains light saturated. However, when the amount of the K d is too high, the total PAR is greatly reduced, and both the diffuse and direct PAR components are reduced (Figure 7). Even though the LUE still increases, the available light is reduced, so the canopy GPP also decreases.

Response of Terrestrial Ecosystem Carbon Fluxes to K d
The terrestrial ecosystem's carbon cycle, especially ecosystem photosynthesis and respiration, is quite sensitive to aerosol loading and cloud fraction (Doughty et al., 2010). Figure 9 shows the relationship among the K d and terrestrial ecosystem carbon fluxes, including GPP, NEE and ER. NEE represents the difference between gross photosynthetic uptake and respiratory losses. The ecosystem is a net carbon sink when the NEE is less than 0, and it becomes a carbon source when the NEE is larger than 0. Therefore, the higher the negative NEE value is, the more carbon is absorbed by the ecosystem. Ecosystem respiration is a biological process that converts organic carbon into CO 2 . It is a compound flux consisting of the aboveground respiration of leaves and woody tissues, the underground respiration of roots (autotrophic soil respiration), and the underground respiration of soil organisms (heterotrophic soil respiration). ER represents a fundamental biospheric function and plays a major role in the global carbon cycle and growth rate of atmospheric CO 2 .
As the aerosol loading and cloud cover increase, the K d increases continuously, and the negative NEE value shows a consistent trend with the GPP; that is, the negative NEE value first increases and then decreases. This shows that the carbon uptake of the ecosystem first increases and then decreases, but when the K d approaches 1, the NEE of the four sites appears to be higher than 0, demonstrating that the ecosystem has switched from a carbon sink to a carbon source. The increase in carbon uptake is primarily a result of the diffuse fertilization effect. The photosynthesis and LUE of shaded leaves are greatly promoted by diffuse PAR. When the K d is too high (around K o ), clouds and aerosols cause a large reduction in total PAR. There is not enough PAR for photosynthesis, so the GPP and NEE decrease. The K d values corresponding to the maximum NEE were 0.56, 0.57, 0.59, and 0.42 at the Cha, Qia, Din and MPM sites, respectively. Except for at the MPM site, where the K d value was lower than K o , the K d values at the other sites were larger than K o , indicating that the NEE response to K d was later than that of the GPP at most sites. This delay effect is caused by ecosystem respiration. The ER response to K d has great uncertainty. The K d has little influence on ER, and ER shows a weak downward trend with increasing K d . This is mainly because the factors influencing ER are complex and diverse, among which the K d is not the leading influencing factor.

The Impact of Environmental Factors on Carbon Flux
The ecosystem carbon cycle consists of biological processes such as photosynthesis and respiration, which are not only affected by PAR and diffuse fertilization effects but also regulated by various environmental factors. Temperature and moisture are important conditions for vegetation growth, so we considered the impact of environmental factors such as temperature (TA), atmospheric vapor pressure deficit (VPD) and latent heat (LE) on the carbon cycle. The result is shown in Figure 10.
Temperature had a significant impact on the GPP of the Cha, Qia, Din, and MPM sites. Among them, the TA of the Qia, Din and MPM sites had an obvious single-peak change in the impact on vegetation, which increased in GPP at first and then decreased. However, the single-peak characteristics at the Cha sites were not significant, showing that the GPP increased with increasing TA. The Cha sites are located in northeastern China, where thermal conditions restrict vegetation growth. Studies have shown that temperature is the main environmental control factor for the increase in GPP in northern Eurasia (Dass et al., 2016). The Qia, Din and MPM sites have lower latitudes and better thermal conditions. Due to higher temperatures there, there was a single-peak change in the impact on vegetation. Studuies have shown that temperature and VPD can limit photosynthesis and respiration (Moazenzadeh et al., 2018). Some researchers believe that vegetation growth has an optimal temperature, and at the most suitable temperature, productivity is the highest. Gui et al. (2021) pointed out that when considering the coeffects of the VPD and TA, the strongest positive effect of diffuse fertilization effects was found at 0-5 hPa and 20-25°C. Zhang Y. et al., 2020 further found that 5-20°C and a VPD < 10 hPa were the optimum ranges for the strongest positive effects of diffuse light on canopy photosynthesis, reflecting that the effect may decrease with increasing stomatal resistance to leaf carbon uptake. Our research shows that when the TA was 26-27°C at the Qia site, the GPP was the highest, the Din site had the highest GPP at approximately 31-32°C, and the MPM site had the highest GPP at approximately 30-31°C. When the temperature is too high, it affects the activity of vegetation enzymes, which in turn affects photosynthesis (Kanniah et al., 2011). At the Cha, Din, and MPM sites, the increase in temperature promoted ER. At the Din and MPM sites, with the increase in temperature, NEE also showed a negative increasing trend, indicating that the increase in temperature at these two sites was conducive to carbon uptake.
VPD is identified as a critical driver of ecosystem function, and the VPD affects the stomatal conductance of the canopy. Similar to the results of previous studies, the GPP we obtained varied with VPD as a unimodal curve, which was mainly manifested at all the FLUXNET forest sites. Studies have shown that when the VPD threshold is 6-8 hPa, and it determines whether the GPP response transitions from a significant increase to a decrease . Our results suggested that an increasing VPD under climate change could reduce forest CO 2 uptake (Sulman et al., 2016). The increase in VPD may exceed the threshold of stomatal openings, which may have a negative impact on GPP. Our research shows that when the VPD of the Cha site was 10-11 hPa, the GPP was the highest, while the Qia and Din sites were 15-16 hPa, and the MPM site was 9-10 hPa. When the VPD is within the appropriate threshold, GPP will increase. However, when the VPD is large, plants will close part of their stomata to reduce the loss caused by transpiration. Studies have shown that in forest ecosystems, GPP is not always restricted by the VPD, which is mainly due to forests having a stronger water holding capacity and deeper root systems (Beer et al., 2007;Chen et al., 2021;Gui et al., 2021). The GPP in these places is mainly affected by radiation or temperature. At the Cha and Din sites, the VPD was helpful in promoting ER, but it was not significant at the other sites. In contrast to the GPP, the impact of the VPD on NEE was U-shaped at the Cha, Din, and MPM sites but it was essentially the same. As the VPD increased, carbon absorption first increased and then decreased. LE is the heat exchange of moisture between the underlying surface and the atmosphere, it has important effects on material and energy fluxes between ecosystems and the atmosphere. High LE can provide energy, water and heat exchange for vegetation (Kasurinen et al., 2014). At all FLUXNET forest sites, the GPP and LE showed a linear increase relationship, while the NEE showed a linear decrease. Among these sites, this pattern was particularly prominent at the Din and MPM sites. These two EBF sites are located in subtropical regions. The vegetation evapotranspiration at these sites is high because of their dense and complex canopy structure. Therefore, the vegetation growth at these sites had a high demand for water and heat. The LE of the MPM site was significantly higher than that of the other sites. The MPM site is dominated by mangroves. Mangroves are tropical wetland vegetation that grow in a hightemperature and humid environment and have high requirements for water and heat, so the latent heat flux deeply influences carbon fluxes (Lugo, 1987). However, the effect of LE on ER is not obvious in general. But at the MPM site, the LE had a positive effect on ER.

CONCLUSION
In this work, the effect of diffuse radiation on the terrestrial ecosystem carbon fluxes of FLUXNET forest sites was investigated. First, we established a K d empirical model based on K t to simulate the K d at the Cha, Qia, Din and MPM sites, and the results showed that the K d model was accurate in estimating the K d (R 2 = 0.88-0.93). The monthly variation in K d was different among the FLUXNET sites. The spatial distribution of K d showed a latitudinal zonality, decreasing gradually from low latitudes to high latitudes. Furthermore, we found changes exist in the relationship between K d , PAR, and GPP and the stage characteristics of GPP. With increasing K d , GPP first increased and then decreased. When K d <K o , the direct PAR decreased rapidly, while the diffuse PAR increased rapidly. When K o <K d <K diff-max , the direct PAR still decreased, and the diffuse PAR still increased, but the total PAR was greatly reduced by excessive aerosols and clouds, so the GPP declined. When K d >K diff-max , the diffuse PAR began to decrease, and the reduction of the superimposed direct PAR caused the GPP of the canopy to drop rapidly. Moreover, under conditions dominated by diffuse light, the LUE of the vegetation canopy was higher than that under direct light. In addition, the results indicated that the ecosystem switches from a carbon sink to a carbon source when the K d approaches 1. Finally, we found that the impact of TA and VPD on GPP is unimodal, and correspondingly, the impact on NEE is U-shaped. However, LE is a combination of hydrothermal conditions that play a positive role in ecosystem productivity and carbon uptake.
Previous studies mostly focused on the effect of the K d or PAR on GPP, and few studies could explain the mechanism by integrating the three factors. Our research clarified the relationship among the K d , PAR, and GPP, as well as K d , LUE, and GPP through mechanistic analysis. It indicated that the diffuse PAR does not always increase with increasing K d but decreases when the K d is too high. In conclusion, the net effect on canopy carbon uptake depends on a balance between the reduction in total PAR and the change in diffuse PAR. However, the K d model based on K t needs to be further improved. In previous studies, it has been concluded that the estimated performance of K d model coupled K t and S is the best. Whereas, due to the limitations of the existing observation data of FLUXNET sites, we failed to establish the K d model of K t and S coupled. In future, we will obtain sunshine duration through modeling or other methods to further improve our model. In addition, the K t that input K d model is also calculated based on the empirical model, because the data used to calculate K t is difficult to be measured due to technical limitations, so the uncertainty may affect the research results (El Mghouchi et al., 2016;Zhang Y. et al., 2017;Shen et al., 2022a). Forests are an important part of the carbon cycle in terrestrial ecosystems. Clarifying the impact of diffuse radiation on the carbon assimilation of forest ecosystems will help improve forest management and promote forest sinks to better achieve the goal of carbon neutrality.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.