Spatially Explicit Modeling of Coupled Water and Carbon Processes Using a Distributed Ecohydrological Model in the Upper Heihe Watershed, China

A fully coupled simulation of ecophysiological, hydrological and biochemical processes is significant for better understanding the individual and interactional impact of sophisticated land surface processes under future disturbances from nature and human beings. In this study, we spatially explicitly modelled evapotranspiration (ET) and photosynthesis (GPP) using a distributed hydrological model, Dynamic Land Model DLM-Ecohydro, over the Upper Heihe watershed for the years of 2013 and 2014. After considering the lateral water movements, the model fairly captured the variations in ET (R2 = 0.82, RMSE = 1.66 mm/day for 2013; R2 = 0.83, RMSE = 1.53 mm/day for 2014) and GPP (R2 = 0.71, RMSE = 5.25 gC/m2/day for 2013; R2 = 0.81, RMSE = 3.38 gC/m2/day for 2014) compared with the measurements from the Arou monitoring station. Vegetation transpiration accounted for total ET of around 65% and 64% in 2013 and 2014, respectively. A large spatial variability was found in these two indicators (14.30–885.36 mm/year for annual ET and 0–2174 gC/m2/day for annual GPP) over the watershed. Soil texture and vegetation functional types were the major factors affecting ET and GPP spatial variability, respectively. The study manifested a coupled water–carbon mechanism through the strong linear relationship between the variations in ET and GPP and the control of hydrological processes on the carbon cycle at the watershed scale. Although the model had a reasonable performance during most parts of the growing seasons, the lack of a soil freezing–thawing scenario caused inevitable discrepancies for the simulation of soil water and heat transfer mechanisms, hence inaccurately estimating the biophysiological processes in the transition period of winter to spring, which should be further improved especially for alpine regions.


Introduction
The land surface ecological processes significantly affect the cycling of energy, mass and moisture between the soil, vegetation and atmosphere. These dynamics also make the land surface an important energy and moisture source, sink and reservoir, and hence contributed a lot to terrestrial and climate systems [1,2]. However, the sophistication of individual, and interactions between, hydrological, ecophysiological and energy cycles are hard to comprehend because of the related nonlinear feedback mechanisms, especially for landscape-scaled ecohydrological processes and the uncertainty of natural and anthropogenic disturbances [3]. The tight linkage between water and carbon dynamics is expressed in many hydrologically controlled ecological processes, for instance, respiration [4], photosynthesis [5], water use efficiency [6] and crop coefficient analyses [7].
Currently, there are many ecophysiological modeling researches published at multiple spatio-temporal resolutions and spans, including point scales [8,9], mesoscales [10,11] and global scales [12] at yearly [13] or minute time-steps [1], while few of them completely considered the soil water processes and the non-linear relationship between it and vegetation ecophysiological processes. Due to the limited understanding of water-carbon interactions, many researches (regardless of modeling or measurement) often regard them as an isolated process [2]. At mesocales, especially for the watershed scale, where the land surface ecological processes are highly sensitive to local water balance, topography is an important contributor for local water cycling and, furthermore, ecophysiological processes. Many researches [14][15][16] suggested that the overlook of topographical-controlled lateral water fluxes would cause significant biases in local estimation of water (evapotranspiration, runoff and soil moisture) and carbon fluxes (GPP: Gross primary productivity, NPP: Net primary productivity and NEP: Net ecosystem productivity). For the traditional hydrological models, although they can completely consider the watershed hydrological processes, they always neglect or simplify the impact of vegetation on soil water and their ecophysiological processes.
Evapotranspiration (ET) and photosynthesis (as GPP) are two critical ecophysiological processes for the connection between land surface water, energy and carbon cycling, and both of them are highly governed by local soil hydrological processes. For GPP, the surrounding water content of the rhizome controls photosynthesis efficiency through the effect of vegetation stomatal and enzyme activity [17]. ET has a stronger interaction with the local water and energy balance. About 70% of the precipitation that falls on the land surface gets back to the atmosphere by vegetation ET, and in arid areas ET loss can be up to 90% [18]. Simultaneously, because vaporization and sublimation would take the absorption of heat, it is also important for local thermal balance. The coupling of water and carbon cycling can better explain these ecohydrological processes and, therefore, help make more efficient decisions for local environmental management.
The Heihe River is the second largest inland river in northwest China and originates from the Qilian Mountains which is in the northern margin of the Tibetan Plateau. It is located in the center of the Eurasian continent, playing an important role in the ancient Silk Road [19]. The upper stream of the watershed takes up around 7% of the Heihe basin area, and is the water source of the whole basin with a wide distribution of grassland which plays a significant role for local carbon assimilation and water cycling [20]. Influenced by the westerly belt circulation and the cold continental polar air mass, the Upper Heihe watershed has an arid to semiarid climate with scarce and concentrated precipitation [21,22]. The water shortage and complex topography are major reasons for sparse vegetation coverage, and hence a more fragile ecological environment [23]. Hitherto, the estimations of ET and GPP in the Upper Heihe watershed were isolated. Yang [24] utilized micro-lysimeters to test the performance of four energy balance methods for ET estimation in the meteorological station of Qilian Mountains, and found that the FAO-56 Penmen-Monteith provided the highest accuracy. Wang [25] established the relationship model between four stations' eddy-covariance-based GPP and the leaf area index (LAI), with vegetation indexes extracted from Moderate Resolution Imaging Spectroradiometer (MODIS) in the Heihe River basin. Yan [20] utilized a remote sensing-based light use efficiency model (MODIS MOD_17) to simulate the GPP of the Upper Heihe watershed, and found that temperature and vapor pressure deficiency had a positive relationships with GPP. These researches lack the description of an ecological water and carbon-coupled mechanism, neither considering the effect of rhizosphere soil water on carbon assimilation rate, nor considering the feedback of the soil carbon pool, biological carbon pool and vegetation growth on the ecosystem water cycle. Sun [26] compared the performance of two land surface models-Community Land Model 4.0 (CLM 4.0) and Dynamic Land Model 1.0 (DLM 1.0)-in the simulation of soil organic carbon of the Upper Heihe watershed, and found that these land surface modes had low accuracies for soil moisture estimation in this cold alpine region. The coupling of the distributed-hydrological and traditional biological models is needed for the watershed scale ecophysiological simulation, because it considers not only the microclimate and topographically controlled soil water movement, but also vegetation-unique ecophysiological characters.
This study utilized a process-based ecohydrological model (DLM-Ecohydro) to simulate two ecophysiological processes, evapotranspiration and photosynthesis, in the cold alpine region. We will test how the model performed after considering topographically controlled lateral water movement and the coupling of the water cycle and carbon assimilation at daily time steps. Through the study of spatiotemporal differentiation of ET and GPP in the Upper Heihe watershed, we will also explore the driving mechanism of ecosystem water and carbon flux variations.

Study Area
The eco-hydrological simulation focused on the Upper Heihe watershed which covers around 10,009 km 2 in the northeastern part of the Tibetan Plateau, China ( Figure 1). The distributions of precipitation and landscapes have strong spatial heterogeneity with a significant topographic relief (from 1600 to 5100 m above sea level). The annual precipitation ranges from 250 mm in lower altitude areas, to 500 mm in higher zones. Melting water from glaciers and snow effectively compensates surface overland and subsurface runoff during spring and early summer [27]. The Upper Heihe watershed is a classical alpine-cold climate with most of the typical landscapes of cold regions, include alpine grassland, alpine meadow, forest, glacier, barren and so on. According to Qin [28], more than 98% of forests fall in the 2000-4000 m a.s.l. area in the watershed. Glaciers lie in places over 4500 m a.s.l. and takes up around 3.4% of the whole watershed area. Grassland and shrubland are mainly distributed in the 1700-4900 m a.s.l. and 2900-5066 m a.s.l. areas, respectively. The annual temperature of the Upper Heihe watershed is lower than 2 • C. Biogeographically, solar and topographic condition are important factors which control the composition of vegetation species. In terms of grassland, Reaumuria soongorica and Sympegma regelu are the main species on sunny slopes, while for shady slopes, it majorly consists of Polygonum viviparum and Carex atrata [29]. Clay is the major soil type within the watershed.

Model Discription
This study utilized a processed-based model, DLM-Ecohydro, which was developed from the Boreal Ecosystem Productivity Simulator (BEPS)-TerrainLab V2.0 [15] and DLM [30]. It modified the modeling period of BEPS-TerrainLab V2.0 from one year to multi-year and localized the vegetation characters for an arid-cold mountain region. Depending on the input datasets, the model is able to resample and run at any spatial resolution. It is carried out using a "four-leave" approach to separate leaf stratifications according to different soil moisture and sunshine states and then spatially upscaled from the leaf to the canopy level. Figure 2 shows a flow chart of the model. Under the input datasets, the model obtains the necessary meteorological and geological conditions, the biophysiological characteristics of local vegetation and the initial ecohydrological status of the watershed. The input energy and precipitation to ecosystem are intercepted according to a hierarchical order at the canopy scale. Rather than a simple bucket model used in Boreal Ecosystem Productivity Simulator (BEPS), DLM-Ecohydro completely considers the effects of lateral fluxes on the above ecophysiological processes by solving a water balance equation. After considering the interactions between ecosystem energy, biochemical and hydrological processes, the model could fairly represent the dynamic states of ecological and hydrological indicators.

Ecophysiological Processes
The DLM-Ecohydro model separates canopy into overstory, understory, moss and soil surface. According to different light conditions and soil moistures states of the rhizosphere, the leaves of each part were divided into sunlit and shade, saturated and unsaturated states. Taking overstory canopy, for example, the conceptualizations of leaf-specific physiological processes is listed as follows: P o,sun,sat = f x R sun , g s,sun,sat , P o,sun,unsat = f x R sun , g s,sun,unsat , P o,shade,sat = f x R sahde , g s,shade,sat , and where P o is overstory canopy transpiration or photosynthesis from: Sunlit leaves with saturated soil water condition (P o,sun,sat ); sunlit leaves with unsaturated soil water conditions (P o,sun,unsat ); shaded leaves with saturated soil water condition (P o,shade,sat ); and shaded leaves with unsaturated water condition (P o,shade,unsat ). Function f x can be either the Penman-Monteith equation for each sub-component of ET (except for sublimation) or the Farquhar model [31] for GPP simulation. R is overstory net radiation for ET and short wave radiation for GPP calculation for sunlit (R sun ) or shaded (R shade ) leaves; g s , overstory stomatal conductance of H 2 O and CO 2 for ET and GPP, respectively, from: Sunlit leaves with saturated soil water condition (g s,sun,sat ); sunlit leaves with unsaturated soil water conditions (g s,sun,unsat ); shaded leaves with saturated soil water condition (g s,shade,sat ); and shaded leaves with unsaturated water condition (g s,shade,unsat ). The up-scaling scheme employed a specific LAI and root fraction as weighting terms, expressed as follows: where P o is overstory canopy ET or GPP, and LAI is the function of the solar zenith angle to divide sunlit (LAI sunlit ) and shaded (LAI shade ) leaves. µ o represented overstory roots in the unsaturated zone and (1 − µ o ) was that in the saturated zone. The fraction of root µ was the function of the root extinction coefficient and water table depth based on the specific pixel and Julian day [32]. The upscaling scheme defaulted the understory canopy just as has shaded leaves. There was no upscaling of soil surface physiological processes because the model deemed it as a single layer with ground specific radiation and soil conductance.

Water Balance Calculation
In the DLM-ecohydro model, the detailed water balance equation considers the lateral water flux by the description of surface overland floor (F surface ) and sub-surface baseflow (F subsurface ) at the pixel scale and daily time steps. Once the infiltration process reaches its max capacity [33], the ponded water would appear on the soil surface and move laterally as F surface . The topography derived lateral subsurface baseflow is the difference between subsurface inflow and outflow within eight neighboring pixels as described in [15]. Water balance calculations are composed with storage changes at unsaturated (∆W unsat ) and saturated (∆W sat ) soil parts which are presented below: where I is total infiltrated water which is comprised of through-fall, snowmelt and stem flow [34], P is percolation and represented the water in the unsaturated zone downward to the saturated part under the impact of gravity, while C is capillary rise upward from the saturated to vadose zone appearing under extreme surface dryness. T unsat and T sat are the transpirations from roots situated in the unsaturated and saturated zones, respectively. E floor is soil evaporation and sublimation of snowpack. The complete soil profile storage change is the sum of two soil zones (Equation (5)).

Spatial Datasets and Model Initialization
Most of the spatial datasets were generated based on the combination of GIS, remote sensing, field monitoring and modeling. All of these datasets were resampled to the same spatial resolution (100 m × 100 m) and projected in World Geodetic System-1984 Coordinate System (WGS-84). The digital elevation model (Figure 3a) was obtained from a 30 m Advanced Spaceborne Thermal Emission and Reflection Radiometer Global Digital Elevation Model (ASTER GDEM) dataset. It was crucial for the computing of lateral hydrological processes and spatial interpolation of various meteorological parameters. Watershed land cover (Figure 3b) was obtained from a 30 m global LC dataset, GlobaLand30, and merged as part of the classifications to fit DLM-Ecohydro. After model parameters localization and optimization, the biophysical parameters for the major land cover types in the Upper Heihe watershed are shown in Table 1 Table 2. The watershed annual maximum LAI (Figure 3d) was generated from Moderate Resolution Imaging Spectroradiometer (MODIS) LAI production (MOD15A2), July 2013. In addition to the spatial datasets above, an initial water table depth (WTD) and ecosystem C-pools (including four biomass C-pools and nine soil C-pools) were necessary for the eco-hydrological and biogeochemical module of DLM-Ecohydro. We used CLM4.0 [35] to spin up WTD based on 2013 meteorological and topographical data, and similar iteratively running processes for current sizes of ecosystem carbon pools were based on the CLM4.0 biogeochemical module. The spun up WTD and ecosystem carbon pools were input as 2013 initial datasets and subsequently used for two years of simulation.

Meteorological Datasets and Flux Measurement
The input meteorological datasets included daily maximum, minimum and average temperature ( • C), median dew point ( • C), total precipitation (mm), wind speed (m/s) and total short-wave radiation (W/m 2 ). These datasets were derived from the Heihe Warershed Allied Telemetry Experimental Research (HiWATER) at the Arou Superstation, which is located in the southeastern part of the study area (100.46 • E, 38.05 • N) [36,37]. Except for the input datasets above, the energy and carbon fluxes and soil moisture data were also prepared for model validation. Simulated volumetric soil water content (VSWC) were compared with a 0-40 cm average measurement with cosmic-ray probes buried at the Arou station, with a footprint scale of 700 m in diameter [38]. The eddy covariance (EC) technique was used in the Arou footprint tower for continuous and long term monitoring of energy and CO 2 fluxes. GPP was derived from the fluxes of net ecosystem exchange (NEE) measured at the EC tower. Based on the Lloyd-Taylar model [39], the relationship between temperature and respiration during the night was used to estimate the daylight respiration rate, thereafter GPP flux could be detected from monitored NEE. The evapotranspiration was converted from measured latent heat fluxes (LE, W/m 2 ) as shown below, where λ is 2.501 × 10 6 J/kg from [35].

Temporal Patterns of Hydrological Indicators
Volumetric soil water content (VSWC) was affected by precipitation, melting runoff, soil temperature and other hydrological or meteorological factors. The model was able to capture the reasonable magnitude and fluctuation of VSWC during most parts of the growing seasons. and lateral water fluxes were the major reasons for soil water losses. There was a dense precipitation after DOY 230 which made the VSWC remain over field capacity. On the other hand, 2013 experienced a prolonged drier autumn that induced a gentle but continuous fall down of VSWC until below the permanent wilting point. While the model overestimated the VSWC during spring and winter when the soil was in a frozen state, it was probably because the model did not take the soil freezing-thawing process into consideration. During the soil frozen period, the permeability would decrease and soil hydraulic characteristics would consequently change. Rather than infiltrated into soil profile, most of the precipitation would drain out as surface overland flow. The lack of these dynamic transformations in the model domain might trigger an unrealistic estimation for soil moisture state. The study area is located in the cold alpine region where the negligence of melting runoff from glaciers and accumulated snowpack might cause the underestimation of soil water especially during spring and early summer.  Figure 5 shows the temporal patterns of two years of lateral water fluxes in the Arou site. On account of the lack of validation data, we only compared it with the daily precipitation. The model captured a rapid response of lateral water fluxes to precipitation in spring and early summer. Under a denser precipitation, the lateral water fluxes of 2014 were much higher than that of 2013. Because in DLM-Ecohydro the surface overland flow (SOLF) was estimated only when saturated water diffuses up to the land surface, i.e., WTD ≤ 0, subsurface baseflow ( SSB) took the major part of the lateral hydrological fluxes. It could reach to 9 mm/day in 2013 and 10 mm/day in 2014 during May and June, respectively. This simulation is more suitable for the soil texture with high hydraulic conductivity, while for clay, it might underestimate the SOLF during dense rainfall periods when the infiltration rate is much lower than the precipitation rate.

Spatial and Temporal Patterns of ET
The temporal patterns and magnitudes of ET are well simulated through the DLM-Ecohydro model, which were able to explain 82% (RMSE = 1.66 mm/day) and 83% (RMSE = 1.53 mm/day) of the variability in 2013 and 2014, respectively. The model agilely captured the seasonal patterns of ET losses at the Arou station of both years: With a fluctuated rise during the early summer until up to the top magnitude (around 5 mm/day) at mid-summer and then decrease to lower than 1 mm/day within September and October (as shown in Figure 6). Although the model could capture most of the variation during a year, there was a general underestimation especially in the winter and early summer. One possible reason for winter ET underestimation was the snow sublimation might cause an increase of latent heat measured through the EC tower. Snowfall is a crucial part of the ecosystem water cycle in the Tibetan Plateau of China [40], which has quite a number of accumulated snowpack during winter. According to the findings of Molotch [41], there could be a rapid response of the EC tower to detect the intercepted snow sublimation from the upper canopy. This contributed to a higher observed ET that the model underestimated. At the same time, the only water input of the DLM-Ecohydro model was precipitation which converts to rainfall or snowpack accumulation depending on air temperature. However, above 4500 m a.s.l. alpine region, ice and snowpack covered all year, which induced an increased snowmelt runoff during spring and summer. Li and Wang [42], utilizing snowmelt runoff model (SRM), demonstrated that May and June are major melting months, and the melting runoff can take more than 50% of the total runoff in the Upper Heihe watershed. The negligence of melting runoff from high altitude snow/ice pack might contribute to the considerable underestimation during early summer. One the other hand, the model calculated plant stomatal conductance with the ambient air temperature, while considering the freezing-thawing process in alpine cold regions, the change of soil temperature usually lagging behind that of air temperature, with a relatively gently periodic variation [43]. It might cause an uncertainty in simulation of stomatal conductance and further cause the discrepancy of ET estimation. For this problem, Govind [34] provided a rectified Jarvis algorithm, including a soil temperature scalar as an environmental factor to calculate stomatal conductance. For further development of the DLM-Ecohydro model, it is necessary to take this modification into consideration. Figure 7 reveals simulated sub-components of Arou station ET in two years. Because the underlying surface of the station is grassland (understory LAI was set as 0), only overstory canopy was took into consideration. In both years, canopy transpiration was the most important component, and took around 65% and 64% of ecosystem ET for 2013 and 2014, respectively. Previous research also demonstrated that ET loss from transpiration can reach 60-80%, especially during the growing season [44]. The flowing loss forms were soil evaporation and canopy intercepted evaporation, accounting for 22-24% and 11-14% in both study years, respectively.  The findings of the simulation showed that the annual ET of the Upper Heihe watershed was around from 14.3 to 885 mm in two study years ( Figure 8). Theoretically, the differences in ecophysiological properties (e.g., g smax , Ω and vegetation height) between vegetation functional types might highly govern the ET losses distribution through their effects of transpirational loss [45]. At the same time, the architectures of different vegetation canopies would affect the over/understory radiation transfer mechanisms [46]. However, spatial heterogeneity of annual total ET had a huge diversity even within the same vegetation stand in the Upper Heihe watershed. The spatial average values of annual ET in deciduous and conifer stands had little discrepancy (131-145 mm) in both years. Annual potential ET losses in the watershed were much more reliant upon soil texture, topography and LAI distribution. One reason for the minimum ET loss of clay was because it had a lower vertical hydraulic conductivity that restricted the water infiltration from soil surface and further limited vegetation transpiration rate. And from HWSD, soil depth of clay was much shallower than that of others. As a result, under the same amount of water input, the rhizosphere of clay was more likely to be located in the saturated zone in the modeling domain. This could cause a lower ET efficiency because the transpirational losses in the saturated zone are much lower than that in unsaturated zone, as shown in Figure 7.
Additionally, annual ET was also highly in concordance with the DEM and LAI distributions in the simulation. Savenije [47] pointed out that topography is one of the most important factors affecting hydrological behavior in watershed scale, and closely related with distribution of plant functional types, soil properties and climate factors. The watershed with a relatively higher ET loss (>200 mm) is mainly concentrated in the south-eastern part of the watershed and near the river outlet area, which has a relatively lower altitude and gentle slope. Comparing with a high elevation mountainous area, the climate of a lower and outlet contiguous area was much warmer and wetter. However, it is important to note that the spatial distribution of soil texture and LAI were obtained from 1 km datasets. Although they were upscaled to 100 m to fit the model, they cannot provide such detailed spatial information. Considering the high dependence of ET on these two variables, the coarse spatial datasets might lead to unrealistic driving forces to some extent. For instance, the inaccurate classification of clay might over/underestimate total ET because of its restricted effect on vegetation transpiration. The model had a good performance in the estimation in GPP magnitudes of two study years. Nevertheless, comparing with measurements, there was a general overestimated GPP in spring and autumn while during June and July there was an underestimation in both years. Li [6] explored the spatiotemporal variations of water use efficiency (WUE), which represents the ratio of CO 2 assimilation to vegetation transpiration [48] in the Heihe River Basin, and pointed that especially in the upper watershed, summer has the highest value of WUE, which is around 1.0-1.2 gC/mm/H 2 O/m 2 comparing with the other seasons (0-0.4 gC/mm/H 2 O/m 2 ). It validated temporal congruency between ET and GPP in the underestimated summer period when the melting runoff was possibly underestimated. From Figure 8b we can see there was a general overestimation of GPP within two years' estimation. This discrepancy might be mainly attributed to the inexact representation of vegetation physiological characteristics consequent on the soil moisture statuses. As discussed above, due to the lack of simulating the soil freezing-thawing process, the model treated each soil texture in the watershed with constant hydraulic properties and remained in the optimal soil moisture condition whole year round. It triggered an overestimated GPP through two major reasons. One is the unrealistic higher stomatal conductance which would speed up photosynthetic process [49]. According to the Jarvis algorithm [50] embedded in DLM-Ecohydro, VSWC has a positive relationship with plant stomatal conductance when it is between the field capacity and permanent wilting point. Actually, for both frozen and thawing period, the infiltrated rate would not be as high as the model representation and quite a number of precipitation or snow melting water probably concentrate on the soil surface and drain out as lateral water fluxes, which would contribute to an increase in soil surface evaporation. This might explain why there was rare overestimation for simulated total ET although it is under the same g s as the GPP calculation. The experiment from Sonnentag [51] also manifested that compared with GPP, there was less sensitivity for ET in different soil moisture conditions. Another possible reason for the bias of modeled GPP was an inaccurate representation of ecosystem nitrogen (N) availability. Via the overestimation of g s , the model exaggerated local photosynthesis, which would raise the soil respiration and soil carbon decomposition and further created a conducive status for increasing biomass and N mineralization. It synergistically promoted photosynthesis and constituted a regional feedback mechanism [52][53][54]. On the other hand, the overestimated soil moisture also positively affected the leaf N content and further overstated GPP through the increased maximum carboxylation rate [55]. That is because deciduous stands had a higher stomatal conductance (g s,max = 5mm/s) and leaf N content (N max = 1.8%) than that of conifers (g s,max = 1.6mm/s; N max = 1.6%). Additionally, conifer canopy was more densely clumped than deciduous stands, with a clumping index of 0.5 and 0.8, respectively. It would decrease the ratio of sunlit leaves in the conifer canopy. Grassland had the most wide distribution within the Upper Heihe and also had the most spatially heterogeneous annual GPP, which ranged from 1.76 to 2119.18 gC m 2 /year in 2013 and 1.75 to 2174.46 gC m 2 /year in 2014. The spatial variability of grassland GPP had a high congruency with LAI, which reveals that the leaf abundance also took a significant contribution of grass photosynthetic efficiency, because it controlled vegetation surface energy balance through the interception of water and solar energy. The distribution of GPP depended on complex interaction effects between topographic conditions, plant function types, vegetation biomass, regional meteorology, soil moisture, temperature and N dynamics in the model. The quantitative relationships and feedbacks between these spatial factors should be further explored. Theoretically, the model could run at any spatial resolution depending on the input dataset availability, while the lack of high-resolution soil types and LAI limited the ability of spatial patterns simulation. To fit the fusion of multi-source datasets, the spatial downscaling of DEM and land cover from 30 m to 100 m would also damage the original data qualities and further decrease the ability of model to explain GPP distribution due to the high correlation between these two factors and photosynthesis.

Co-Drive Mechanism of Ecosystem ET and GPP
According to the above research we found that there is a strong synchronism between ecosystem evapotranspiration and photosynthesis. Figure 11 illustrates that the simulated ET and GPP have a significant positive liner relation (R 2 = 0.88) within two study years. Except for a consistent seasonal variation with a same peak period at July and August, they also have a similar spatial heterogeneity. Previous research also pointed the consistent diurnal variation between ecosystem water and carbon fluxes [56], and Steduto [57] found that the annual accumulated carbon assimilation and evapotranspiration have a stable linear relationship. Except for the impact of soil moisture, another major reason is probably because the driving energies of ecosystem C-water cycle have a synchronous change. Canopy intercepted net radiation for ET and short-wave radiation for photosynthesis have a similar co-linear relation with input solar radiation. Additionally, as a common organ for transpiration and photosynthesis, the growth of leaves synchronously controls canopy ET and GPP fluctuation and this control is hardly affected by soil water and fertilizer conditions [58]. It also explains the similar spatial distribution between C-water fluxes and regional LAI discussed above. Furthermore, daily temperature and wind speed can also synchronously affect the exchange rate of CO 2 and H 2 O between vegetation and the atmosphere through the control of enzymatic activity, stomatal conductance and vaporizing rate.

Conclusions
Based on two fundamental ecophysiological indicators, ET and GPP, this study utilized the DLM-Ecohydro model to simulate the exchanges between energy, water and carbon assimilation in the Upper Heihe watershed with a tightly coupled manner. The combination of in situ measurements and remotely sensed datasets provided detailed environmental disturbances from both nature and human beings, contributing a regionally accurate and spatially explicit simulation. The land surface model fully took topographical influence into consideration, modeling the lateral hydrological processes by solving a water balance equation. In order to accurately represent the ecophysiological processes in the watershed, the model separately simulated sunlit and shaded leaves and flexibly derived soil profiles as multiple layers according to the moisture regimes.
The simulation suggested that in spatial and temporal patterns, ET and GPP have significantly synchronous variations in the Upper Heihe watershed. Although the precipitation was mainly concentrated during May and September in the watershed, soil water content was still a major limiting factor for vegetation transpiration and photosynthesis efficiency in mid-summer. The spatial distribution of ET highly depends on soil texture. Silt loam and clay had the maximum and minimum annual ET loss, respectively, mainly because of their different abilities of hydraulic conductance. Vegetation transpiration accounted for the most part of ecosystem ET losses in the two study years. Photosynthesis was more sensitive to plant functional types: Deciduous forests had the highest annual carbon assimilation rate owing to it having a higher stomatal conductance and leaf nitrogen content and a lower canopy clumping index. LAI and topography, however, made a strong contribution for both ecophysiological processes simultaneously.
There are also some limitations of this study. Because the model did not consider the dynamics of frozen soil and the melting runoff from glaciers in the alpine regions, it would cause an unreality estimation of soil moisture regimes during the frozen and major melting seasons. While overall, the model showed a fair performance for the simulation of hydro-ecological processes, especially during vegetation growing seasons.