A One-Source Approach for Estimating Land Surface Heat Fluxes Using Remotely Sensed Land Surface Temperature

Abstract

The partitioning of available energy between sensible heat and latent heat is important for precise water resources planning and management in the context of global climate change. Land surface temperature (LST) is a key variable in energy balance process and remotely sensed LST is widely used for estimating surface heat fluxes at regional scale. However, the inequality between LST and aerodynamic surface temperature (T_aero) poses a great challenge for regional heat fluxes estimation in one-source energy balance models. To address this issue, we proposed a One-Source Model for Land (OSML) to estimate regional surface heat fluxes without requirements for empirical extra resistance, roughness parameterization and wind velocity. The proposed OSML employs both conceptual VFC/LST trapezoid model and the electrical analog formula of sensible heat flux (H) to analytically estimate the radiometric-convective resistance (r_ae) via a quartic equation. To evaluate the performance of OSML, the model was applied to the Soil Moisture-Atmosphere Coupling Experiment (SMACEX) in United States and the Multi-Scale Observation Experiment on Evapotranspiration (MUSOEXE) in China, using remotely sensed retrievals as auxiliary data sets at regional scale. Validated against tower-based surface fluxes observations, the root mean square deviation (RMSD) of H and latent heat flux (LE) from OSML are 34.5 W/m² and 46.5 W/m² at SMACEX site and 50.1 W/m² and 67.0 W/m² at MUSOEXE site. The performance of OSML is very comparable to other published studies. In addition, the proposed OSML model demonstrates similar skills of predicting surface heat fluxes in comparison to SEBS (Surface Energy Balance System). Since OSML does not require specification of aerodynamic surface characteristics, roughness parameterization and meteorological conditions with high spatial variation such as wind speed, this proposed method shows high potential for routinely acquisition of latent heat flux estimation over heterogeneous areas.

1. Introduction

Land surface heat fluxes are essential components of water and energy cycles interactions in the hydrosphere, atmosphere and biosphere [1,2,3,4,5,6]. The energy balance and its partitioning between sensible heat flux (H) and latent heat flux (LE) are extremely important for understanding global climate change and land-atmosphere interaction [7]. Satellite remote sensing provides a promising way to obtain information on both LE and H fluxes at larger spatiotemporal scales than models depend on in-situ meteorological observations alone [8,9,10]. Considerable progress has been made in the development of global flux products from space observations [11,12,13,14,15,16].

Remotely sensed land surface temperature (LST) is widely used in evapotranspiration (ET) modeling, and it is a key variable for partitioning available energy into H and LE [9,17]. LST is often referred to as “radiative surface temperature” [18]. Because of the inequality between LST and aerodynamic surface temperature (T_aero), one-source and two-source surface energy balance models have been proposed with various expressions for aerodynamic resistance, in order to estimate the sensible heat loss based on difference between surface and air temperature [19,20,21,22,23]. Much effort has been made to compensate for the difference between LST and T_aero in one-source models, most of which was concentrated on making corrections by: (1) adjusting temperature based on empirical relationship between LST and T_aero [19,21,23]; and (2) adjusting the aerodynamic resistance roughness length for heat or the kB⁻¹ parameter, or by including empirical “extra resistance” (r_ex) factors [24,25]. For instance, Chehbouni et al. [20] and Matsushima [22] proposed methods to interpret the relationship between LST and T_aero, and make further correction on LST based on this relationship. However, due to the high empirical dependency of such models and the non-uniqueness in the relationship between LST and T_aero, these approaches generally produced significant errors when applied to partially-vegetated landscapes [26,27]. Other works were based on adjustments of the roughness length for heat or the kB⁻¹ parameter, or by including empirical “extra resistance” factors [23,25,28,29]. However, the high empirical dependency of the model parameterization on aerodynamic characteristics introduces great uncertainty in H estimation [9,10,17,24]. Extensive reviews on the methodologies of compensation for the difference between LST and T_aero can be found in the related works [17,22,24,30]. Two-source models, which eliminate the need for extra resistance and treat H as the sum of sensible heat fluxes from both soil surface and vegetation, were developed by Norman et al. [26]. However, various aerodynamic resistances need to be parameterized in two-source models, and the soil surface aerodynamic resistance to heat transfer is the most difficult parameter to be characterized. Extensive reviews of satellite-based methodologies on surface heat fluxes estimation can be found in the related works [10,17,30].

Lacking information on regional wind speed and aerodynamic surface characteristics currently poses a major challenge for regional surface heat fluxes estimation [17,31]. In particular, wind velocity at regional scale is highly variable and is difficult to estimate from in-situ observations [32,33]. In addition, estimation of the aerodynamic resistance usually requires stability corrections, as well as determination of parameters including suitable roughness lengths and zero displacement level [25,34,35,36]. Boegh et al. [31] pointed out that aerodynamic resistance is usually estimated from local data, although it must be noted that areal average of roughness length from point-scale observation is highly unreliable due to its nonlinearity. Several studies have been carried out in an attempt to estimate surface energy fluxes without having to calculate various aerodynamic resistances. Such attempt includes the proposal of the spatial contextual models based on the LST and vegetation fractional cover (VFC) space (hereinafter referred to as VFC/LST space), and a key advantage of this category of approaches is their relative independence from site-specific tuning of resistances-related parameters [36,37,38,39,40,41]. However, a large number of pixels with a wide range of soil wetness and fractional cover are required in such models. They also involve subjectivity in the determination of so-called “dry edge” and “wet edge” [42]. The alternative approach is to solve LE without the need for estimation of the aerodynamic resistance. For instance, Mallick et al. [43] proposed a new method, named Surface Temperature Initial Closure (STIC), which integrates the LST into the Penman–Monteith equation for estimating the terrestrial surface energy balance fluxes. Boegh et al. [31] proposed a novel method that eliminates the need for introducing the parameter of excess resistance and used three equations to estimate atmospheric resistance and surface resistance. However, this method needs a modification of the surface humidity parameter (i.e., h_s,max in Boegh et al. [31]).

The objective of this study is twofold: (1) to develop a One-Source Model for Land (OSML) for estimating land surface fluxes which circumvents the difficulties of parameterizing extra resistance and acquiring regional-scale wind velocity; and (2) to evaluate the performance of the proposed OSML model by validation with in-situ surface fluxes observations and comparison with SEBS in terms of estimating both surface fluxes and radiometric-convective resistance. Section 2.1 and Section 2.2 present a theoretical background and detailed description of the OSML model, and Section 2.3 introduces the commonly used one-source model (SEBS (Surface Energy Balance System)), which served as a benchmark for evaluating OSML in this study. Section 3 introduces the remotely-sensed data sets used for driving OSML and the eddy covariance system measurements used for assessing the model. Section 4 reports the LE and H fluxes validation results from OSML at two experimental regions with contrasting climatic features and different surface complexities, namely, Soil Moisture-Atmosphere Coupling Experiment (SMACEX) site and Multi-Scale Observation Experiment on Evapotranspiration (MUSOEXE) site. OSML performance is also compared with SEBS in estimated surface fluxes. Section 5 discusses the advantages and limitations of this work. Section 6 provides a conclusion.

2. Model Description

2.1. Theoretical Background of One-Source Model

The energy balance equation can be formulated for the whole canopy-soil system:

$${\mathrm{R}}_{\mathrm{n}}=\mathrm{H}+\mathrm{LE}+\mathrm{G}$$

(1)

where R_n is the net radiation flux (W/m²), H is the sensible heat flux (W/m²), LE is the latent heat flux (W/m²) and G is the soil heat flux (W/m²), which can be estimated using remote sensing retrievals [25,44]. The key issue of one-source model is to accurately estimate H, and then LE is calculated as the residual term in the energy balance equation. In one-source model, H is generally calculated using the electrical analog:

$$\mathrm{H}={\mathsf{\rho }\mathrm{C}}_{\mathrm{p}}\frac{{\mathrm{T}}_{\mathrm{aero}}-{\mathrm{T}}_{\mathrm{a}}}{{\mathrm{r}}_{\mathrm{a}}}={\mathsf{\rho }\mathrm{C}}_{\mathrm{p}}\frac{{\mathrm{T}}_{\mathrm{s}}-{\mathrm{T}}_{\mathrm{a}}}{{\mathrm{r}}_{\mathrm{a}}+{\mathrm{r}}_{\mathrm{ex}}}={\mathsf{\rho }\mathrm{C}}_{\mathrm{p}}\frac{{\mathrm{T}}_{\mathrm{s}}-{\mathrm{T}}_{\mathrm{a}}}{{\mathrm{r}}_{\mathrm{ae}}}$$

(2)

where ρ is the air density (kg/m³), and C_p is the specific heat of air at constant pressure (J/kg·°C). T_a is the air temperature (°C). T_aero is the aerodynamic temperature (°C) which is rather difficult to estimate over heterogeneous surface and cannot be directly obtained from remote sensing retrievals. The variable r_a is aerodynamic resistance (s/m).

Remotely sensed land surface temperature (T_s, °C), on the contrary, is more easily obtained and is used as a replacement for T_aero. The difference between T_s and T_aero may differ by several degrees. In order to mediate the difference between T_s and T_aero, an extra resistance (r_ex) is added in the denominator (Equation (2)). The sum of r_a and r_ex is termed as radiometric-convective resistance (r_ae, s/m) by Norman et al. [26]. Boegh et al. [31] termed r_ae as the “atmospheric resistance” [31]. In one-source model r_ex is often calculated using the kB⁻¹ approach [25,29]. However, kB⁻¹ is often empirically determined by local calibration and the erratic values of kB⁻¹ for sparse vegetation can lead to large uncertainties of one-source model when applied to heterogeneous surface. Full reviews of the kB⁻¹ approaches on one-source model can be found in Paul et al. [24] and Boegh et al. [31].

2.2. Development of One-Source Model for Land (OSML)

In this study, the parameterization scheme of r_ae is inspired by the work of Moran et al. [45] and Boegh et al. [31]. Moran et al. [45] introduced a water deficit index (WDI) based on the VFC/LST trapezoid space and extended the application of the crop water stress index (CWSI) from fully- to partially-vegetated surface areas. Based on the WDI, H can be estimated using the following approach:

$$\mathrm{H}=({\mathrm{R}}_{\mathrm{n}}-\mathrm{G})-\mathrm{LE}=({\mathrm{R}}_{\mathrm{n}}-\mathrm{G})-(1-\mathrm{WDI}){\mathrm{E}}_{\mathrm{p}}$$

(3)

$$\mathrm{WDI}=\frac{\mathrm{a}}{\mathrm{a}+\mathrm{b}}$$

(4)

$$\mathrm{a}=({\mathrm{T}}_{\mathrm{s}}-{\mathrm{T}}_{\mathrm{a}})-{\mathrm{f}}_{\mathrm{c}}({{\mathrm{T}}_{\mathrm{s}}-{\mathrm{T}}_{\mathrm{a}})}_{\mathrm{C}}+({\mathrm{f}}_{\mathrm{c}}-1)({\mathrm{T}}_{\mathrm{s}}-{\mathrm{T}}_{\mathrm{a}}{)}_{\mathrm{D}}$$

(5)

$$\mathrm{b}=(1-{\mathrm{f}}_{\mathrm{c}})({\mathrm{T}}_{\mathrm{s}}-{\mathrm{T}}_{\mathrm{a}}{)}_{\mathrm{A}}+{\mathrm{f}}_{\mathrm{c}}({{\mathrm{T}}_{\mathrm{c}}-{\mathrm{T}}_{\mathrm{a}})}_{\mathrm{B}}-({\mathrm{T}}_{\mathrm{s}}-{\mathrm{T}}_{\mathrm{a}})$$

(6)

where a is the temperature difference between Point M (Figure 1) and the cold edge, and b is the temperature difference between Point M and the dry edge. E_p is the potential evaporation and is computed using the Priestley–Taylor equation:

$${\mathrm{E}}_{\mathrm{p}}=1.26\frac{\mathsf{\Delta }}{\mathsf{\Delta }+\mathsf{\gamma }}({\mathrm{R}}_{\mathrm{n}}-\mathrm{G})$$

(7)

where Δ is the slope of the saturation vapor pressure curve at air temperature T_a (kPa·°C⁻¹), and γ is the psychrometric constant (kPa·°C⁻¹). For a given pixel, the four theoretical extreme points can be derived based on the values of R_n, G, T_a, vapor pressure deficit (VPD) and r_ae following the equations from Moran et al. [45]. That is, for the dry bare soil of Point A, evaporation equals zero and the difference between LST and T_a can be derived as follows:

$${({\mathrm{T}}_{\mathrm{s}}-{\mathrm{T}}_{\mathrm{a}})}_{\mathrm{A}}={\mathrm{r}}_{\mathrm{ae}}{(\mathrm{R}}_{\mathrm{n}}{-\mathrm{G})/\mathsf{\rho }}_{\mathrm{a}}{\mathrm{C}}_{\mathrm{p}}$$

(8)

For the saturated bare soil of Point D where the surface resistance equals zero, the difference takes the form of:

$${{(\mathrm{T}}_{\mathrm{s}}{-\mathrm{T}}_{\mathrm{a}})}_{\mathrm{D}}={[(\mathrm{r}}_{\mathrm{ae}}{(\mathrm{R}}_{\mathrm{n}}{-\mathrm{G}))/(\mathsf{\rho }}_{\mathrm{a}}{\mathrm{C}}_{\mathrm{p}})][\mathsf{\gamma }/(\mathsf{\Delta }+\mathsf{\gamma })]-[\mathrm{VPD}/(\mathsf{\Delta }+\mathsf{\gamma })]$$

(9)

For the well-watered full-cover vegetation of Point C,

$$\begin{array}{l}{({\mathrm{T}}_{\mathrm{s}}-{\mathrm{T}}_{\mathrm{a}})}_{\mathrm{C}}=[({\mathrm{r}}_{\mathrm{ae}}({\mathrm{R}}_{\mathrm{n}}-\mathrm{G}))/({\mathsf{\rho }}_{\mathrm{a}}{\mathrm{C}}_{\mathrm{p}})][\mathsf{\gamma }(1+{\mathrm{r}}_{\mathrm{cp}}/{\mathrm{r}}_{\mathrm{ae}})/(\mathsf{\Delta }+\mathsf{\gamma }(1+{\mathrm{r}}_{\mathrm{cp}}/{\mathrm{r}}_{\mathrm{ae}}))]\\ -[\mathrm{VPD}/(\mathsf{\Delta }+\mathsf{\gamma }(1+{\mathrm{r}}_{\mathrm{cp}}/{\mathrm{r}}_{\mathrm{ae}}))]\end{array}$$

(10)

For the water-stressed full-cover vegetation of Point B,

$$\begin{array}{l}{({\mathrm{T}}_{\mathrm{c}}-{\mathrm{T}}_{\mathrm{a}})}_{\mathrm{B}}=[({\mathrm{r}}_{\mathrm{ae}}({\mathrm{R}}_{\mathrm{n}}-\mathrm{G}))/({\mathsf{\rho }}_{\mathrm{a}}{\mathrm{C}}_{\mathrm{p}})][\mathsf{\gamma }(1+{\mathrm{r}}_{\mathrm{cx}}/{\mathrm{r}}_{\mathrm{ae}})/(\mathsf{\Delta }+\mathsf{\gamma }(1+{\mathrm{r}}_{\mathrm{cx}}/{\mathrm{r}}_{\mathrm{ae}}))]\\ -[\mathrm{VPD}/(\mathsf{\Delta }+\mathsf{\gamma }(1+{\mathrm{r}}_{\mathrm{cx}}/{\mathrm{r}}_{\mathrm{ae}}))]\end{array}$$

(11)

where r_cp is the canopy resistance for full cover well-watered vegetation (s/m), and r_cx is the canopy resistance for full cover vegetation with no available water (s/m). Values of r_cp and r_cx could be parameterized following Jarvis-type representation of canopy resistance as follows [46,47,48,49]:

$${\mathrm{r}}_{\mathrm{cp}}={\mathrm{r}}_{\mathrm{cmin}}{\mathrm{F}}_1^{-1}{\mathrm{F}}_2^{-1}{\mathrm{F}}_3^{-1}{\mathrm{F}}_{cp}^{-1}/\mathrm{LAI}$$

(12)

$${\mathrm{r}}_{\mathrm{cx}}={\mathrm{r}}_{\mathrm{cmin}}{\mathrm{F}}_1^{-1}{\mathrm{F}}_2^{-1}{\mathrm{F}}_3^{-1}{\mathrm{F}}_{cx}^{-1}/\mathrm{LAI}$$

(13)

where ${\mathrm{r}}_{\mathrm{cmin}}$ and ${\mathrm{r}}_{\mathrm{cmax}}$ is the minimum value of the canopy resistance. LAI is leaf area index. The factor F₁ describes the influence of the photosynthetically active radiation on canopy resistance; F₂ describes the influence of vapor pressure deficit on canopy resistance and F₃ describe the influence of air temperature on canopy resistance. The factors F_cp and F_cx take into account the effect of water stress on canopy resistance. For canopy resistance for full cover well-watered vegetation, F_cp takes the value of 1.0. For canopy resistance for full cover vegetation with no available water, F_cx takes the value of 0.05. Estimation schemes for F₁, F₂ and F₃ are taken from Dickinson et al. [46], Noihan and Planto [48] and Wigmosta et al. [50] as:

$${\mathrm{F}}_1=\frac{1+f}{f+({\mathrm{r}}_{\mathrm{cmin}}/{\mathrm{r}}_{\mathrm{cmax}})}$$

(14)

$$f=0.55\frac{R_S^{\downarrow }}{R_{SL}^{\downarrow }}\frac{2.0}{LAI}$$

(15)

$${\mathrm{F}}_2=1.0-VPD/c_{vpd}$$

(16)

$${\mathrm{F}}_3=0.08{\mathrm{T}}_{\mathrm{a}}-0.0016{\mathrm{T}}_a^2$$

(17)

where ${\mathrm{r}}_{\mathrm{cmax}}$ is the maximum value of the canopy resistance, $R_S^{\downarrow }$ is the downward shortwave radiation, and $R_{SL}^{\downarrow }$is a limit value of 100 W/m² for crops. $c_d$ is the vapor pressure deficit causing stomatal closure, $c_{\mathrm{vp}d}$ is about 4 kPa.

Equations (2)–(11) are combined to solve r_ae (please see Appendix A for details). The computations are performed using pixel-based LST, R_n, G, T_a, and VPD. The estimated r_ae are used to compute H, and consequently, LE can be solved as the residual term in the energy balance equation (Equation (1)). The model inputs (from meteorological observations and remote sensing retrievals), intermediate variables and the main procedures for implementation of OSML are summarized in a flow chart (Figure 2).

2.3. Surface Energy Balance System (SEBS)

The SEBS model was developed for estimating regional surface heat fluxes by using satellite remote sensing data sets [25]. In this study, SEBS serves as a benchmark for evaluating OSML performance. In SEBS, H is calculated iteratively by solving three nonlinear equations:

$$\mathrm{u}=\frac{{\mathrm{u}}_{*}}{\mathrm{k}}[\ln (\frac{{\mathrm{z}-\mathrm{d}}_0}{{\mathrm{z}}_{\mathrm{om}}})-{\mathsf{\psi }}_{\mathrm{m}}(\frac{{\mathrm{z}-\mathrm{d}}_0}{\mathrm{L}})+{\mathsf{\psi }}_{\mathrm{m}}(\frac{{\mathrm{z}}_{\mathrm{om}}}{\mathrm{L}})]$$

(18)

$${\mathsf{\theta }}_0-{\mathsf{\theta }}_{\mathrm{a}}=\frac{\mathrm{H}}{{\mathrm{ku}}_{*}{\mathsf{\rho }\mathrm{C}}_{\mathrm{p}}}[\ln (\frac{{\mathrm{z}-\mathrm{d}}_0}{{\mathrm{z}}_{\mathrm{oh}}})-{\mathsf{\psi }}_{\mathrm{h}}(\frac{{\mathrm{z}-\mathrm{d}}_0}{\mathrm{L}})+{\mathsf{\psi }}_{\mathrm{h}}(\frac{{\mathrm{z}}_{\mathrm{oh}}}{\mathrm{L}})]$$

(19)

$$\mathrm{L}=-\frac{{\mathsf{\rho }\mathrm{C}}_{\mathrm{p}}{\mathrm{u}}_{*}^3{\mathsf{\theta }}_{\mathrm{v}}}{\mathrm{kgH}}$$

(20)

where u is the wind speed (m/s), ${\mathrm{u}}_{*}$ is the friction velocity (m/s), k is the von Kaman’s constant, d₀ is the zero plane displacement (m), z is the height above the surface (m), and z_om and z_oh are the roughness heights for momentum and heat transfer (m), respectively. θ₀ and θ_a are the potential temperatures at surface and at height z (K), respectively. ${\mathsf{\psi }}_{\mathrm{m}}$ and ${\mathsf{\psi }}_{\mathrm{h}}$ are the stability correction functions for sensible heat and momentum transfer, respectively. L is the Obukhov length (m), and g is the acceleration due to gravity (m/s²).

In SEBS, the roughness height for heat transfer is calculated using a model for calculating kB⁻¹ [25,51]. A weighted average between limiting cases of full canopy, bare soil conditions and mixed vegetation is implemented based on vegetation fractional cover [28,51]:

$${\mathrm{z}}_{\mathrm{oh}}=\frac{{\mathrm{z}}_{\mathrm{om}}}{{\exp (\mathrm{kB}}_{\mathrm{SEBS}}^{-1})}$$

(21)

$${\mathrm{kB}}_{\mathrm{SEBS}}^{-1}=\frac{{\mathrm{kC}}_{\mathrm{d}}}{{4\mathrm{C}}_{\mathrm{t}}\frac{{\mathrm{u}}_{*}}{\mathrm{u}(\mathrm{h})}{(1-\mathrm{e}}^{-\mathrm{n}/2})}{\mathrm{f}}_{\mathrm{c}}^2+\frac{\mathrm{k}\frac{{\mathrm{u}}_{*}}{\mathrm{u}(\mathrm{h})}\frac{{\mathrm{z}}_{\mathrm{om}}}{\mathrm{h}}}{{\mathrm{C}}_{\mathrm{t}}^{*}}{\mathrm{f}}_{\mathrm{c}}^2{\mathrm{f}}_{\mathrm{s}}^2+{\mathrm{kB}}_{\mathrm{s}}^{-1}{\mathrm{f}}_{\mathrm{s}}^2$$

(22)

where C_d is the drag coefficient of the foliage, C_t is the heat transfer coefficient, u(h) is the horizontal wind speed at the canopy height, n is the wind speed profile extinction within the canopy, f_c is vegetation fractional cover, f_s is the fractional cover for soil, ${\mathrm{C}}_{\mathrm{t}}^{*}$ is the heat transfer coefficient of the soil, and ${\mathrm{kB}}_{\mathrm{s}}^{-1}$ is the value for bare soil surface, calculated as:

$${\mathrm{kB}}_{\mathrm{s}}^{-1}=2.46{(\mathrm{R}{\mathrm{e}}_{*})}^{1/4}-\ln (7.4)$$

(23)

where $\mathrm{R}{\mathrm{e}}_{*}$ is the roughness Reynolds number. In this study, we briefly introduced the SEBS algorithm for H estimation, and the detailed description of the model can be found in Su [51].

SEBS requires three sets of input data. The first set is satellite-based land surface variables, including albedo, emissivity, LST, and LAI. The second set includes meteorological variables, such as air temperature, water vapor pressure and wind speed. The third set is atmospheric radiation fluxes, including downward shortwave radiation and downward longwave radiation. In this study, the SEBS module in ILWIS (Integrated Land and Water Information System) was adopted. In order to compare the differences between SEBS and OSML, the radiometric-convective resistance of SEBS (r_{ae_SEBS}) is calculated as:

$${\mathrm{r}}_{\mathrm{ae}\text{\_}\mathrm{SEBS}}={\mathsf{\rho }\mathrm{C}}_{\mathrm{p}}\frac{{\mathrm{T}}_{\mathrm{s}}{-\mathrm{T}}_{\mathrm{a}}}{{\mathrm{H}}_{\mathrm{SEBS}}}$$

(24)

where H_SEBS is the sensible heat flux derived from the SEBS module.

3. Study Area and Data Processing

3.1. SMACEX Campaign

3.1.1. Description of SMACEX Campaign Area and Related Experiment SMACEX

The Soil Moisture-Atmosphere Coupling Experiment (SMACEX) was conducted in conjunction with the Soil Moisture Experiment 2002 (SMEX02) during June and July 2002 near Ames, Iowa. The SMACEX was designed to provide direct measurement/remote-sensing/modeling approaches for understanding the impact of spatial and temporal variability in vegetation cover, soil moisture and other land surface states on the turbulent flux exchange with the atmosphere. A full description of the SMACEX experiment can be found in Kustas et al. [52]. The land cover in the experimental area is primarily corn and soybean. The extensive measurements of atmospheric states, vegetation conditions, soil properties and surface fluxes allow for rigorous analysis and the validation of land surface model and remote sensing observations in the study area [53,54].

3.1.2. Remotely Sensed and In-Situ Observations Required for OSML Input over SMACEX

The required in-situ meteorological forcing and satellite-based retrievals for OSML and SEBS are listed in

Table 1. Meteorological forcing and satellite-based observations required for OSML and SEBS for SMACEX site.

	OSML		SEBS
	DOY174	DOY182	DOY174	DOY182
Meteorological Forcing
Incoming solar radiation (W/m2)	834	859	834	859
Air temperature (°C)	29.6	29.35	29.6	29.35
Vapor Pressure (kPa)	1.96	2.2	1.96	2.2
Atmospheric pressure (kPa)	98.2	98.2	98.2	98.2
Wind speed (m/s)	×	×	6.3	5.32
Satellite-Based Retrievals
LST	Derived from TIR, band 6 based on Li et al. [55]
LAI	Empirical relationship given by Anderson et al. [60]
NDVI	Calculated with NIR and RED band
Land cover	Obtained from the National Snow and Ice Data Center (http://nsidc.org/index.html)
Vegetation Fractional Cover	Estimated from NDVI
Albedo	Albedo was retrieved from the visible and near-infrared bands (1–5, 7) of the Landsat images following Tasumi et al. [57]
Crop height (m)	×	×	Empirical relationship given by Anderson et al. [60]
Surface roughness (m)	×	×	0.1 × Crop height
Displacement height (m)	×	×	0.67 × Crop height

× means not required.

. Two Landsat scenes were used: one from the Landsat-5 Thematic Mapper (TM) on 23 June 2002 (Day of Year (DOY) 174), and the other was on 1 July 2002 (DOY 182) from the Enhanced Thematic Mapper plus (ETM+) on Landsat 7. The two Landsat scenes were obtained from the National Snow and Ice Data Center (http://nsidc.org/index.html).

The LST was derived from the thermal infrared band (TIR, band 6). The atmospheric transmittance, spectral radiance and downwelling sky radiance from the atmosphere were provided by Li et al. [55]. The visible, near-infrared (NIR) bands were re-sampled to be consistent in resolution with images from the thermal infrared band. Details concerning the atmospheric correction and retrieval of LST were given in Li et al. [55] and Li et al. [56]. Albedo was retrieved from the visible and near-infrared bands (1–5, 7) of the Landsat images following Tasumi et al. [57] and Allen et al. [58,59]. The land cover classification for the study area was obtained from the National Snow and Ice Data Center (http://nsidc.org/index.html). The vegetation fractional cover was estimated from the normalized difference vegetation index (NDVI) following:

$${\mathrm{f}}_{\mathrm{c}}=1-{(\frac{{\mathrm{NDVI}-\mathrm{NDVI}}_{\min }}{{\mathrm{NDVI}}_{\max }{-\mathrm{NDVI}}_{\min }})}^{\mathrm{k}}$$

(25)

where NDVI_max and NDVI_min are NDVI values for complete vegetation cover and bare soil, respectively. The coefficient k is a function of the leaf orientation distribution within the canopy, which is set as 0.6175 in this study, as suggested by Li et al. [56].

In order to fairly compare the performance of OSML and SEBS, the same net radiation (R_n) and soil heat flux (G) were used for the two models. The net radiation R_n is estimated as:

$${\mathrm{R}}_{\mathrm{n}}={(1-\mathsf{\alpha })\mathrm{S}}_{\mathrm{wd}}+{\mathsf{\epsilon }\mathsf{\delta }(\mathsf{\epsilon }}_{\mathrm{a}}{{(\mathrm{T}}_{\mathrm{a}}+273.15)}^4-{{(\mathrm{T}}_{\mathrm{s}}+273.15)}^4)$$

(26)

where α is the albedo, ε is the emissivity, S_wd is the downwelling shortwave radiation (W/m²), and δ is the Stefen–Boltzman constant. The soil heat flux G is parameterized following Choudhury et al. [44]:

$$\mathrm{G}={\mathrm{R}}_{\mathrm{n}}{(\mathrm{f}}_{\mathrm{c}}{\mathsf{\tau }}_{\mathrm{v}}{+\mathsf{\tau }}_{\mathrm{s}}(1-{\mathrm{f}}_{\mathrm{c}}))$$

(27)

where τ_v and τ_s are the soil heat flux ratios for full vegetation (0.05) and bare soil (0.315).

Empirical relationships between the Landsat-derived vegetation index and LAI, vegetation cover fraction and crop height over the watershed were given by Anderson et al. [60]. In this study, these empirical relationships were used to derive LAI and vegetation height for the SMACEX region. For both corn and soybean, LAI was derived using the following equation:

$$\mathrm{LAI}=(2.88\mathrm{NDWI}+1.14)(1+0.104\exp (4.1\mathrm{NDWI}))$$

(28)

where NDWI is the normalized difference water index, which can be computed from the near-infrared (NIR, band 4) and shortwave infrared (SWIR, band 5) reflectances:

$$\mathrm{NDWI}=(\mathrm{NIR}-\mathrm{SWIR})/(\mathrm{NIR}+\mathrm{SWIR})$$

(29)

Crop heights for corn (h_corn) and soybean (h_soybean) were derived using the following equations:

$${\mathrm{h}}_{\mathrm{corn}}=(1.2\mathrm{NDWI}+0.6)(1+0.04\exp (5.3\mathrm{NDWI}))$$

(30)

$${\mathrm{h}}_{\mathrm{soybean}}=(0.5\mathrm{NDWI}+0.26)(1+0.005\exp (4.5\mathrm{NDWI}))$$

(31)

3.1.3. Eddy Covariance System Measurements over SMACEX

Twelve towers with eddy covariance (EC) systems were deployed in the central area of the experiment. All towers were instrumented with sensors to measure turbulent fluxes and soil heat fluxes, six of them are corn sites and six soybean sites. The land use types of study area and the spatial distribution of flux towers are shown in Figure 3. It is well recognized that the field measurements of heat fluxes using the eddy covariance technique often fail to show closure of the surface energy budget [61,62,63]. During SMACEX, the energy balance closure was 0.71 for the soybean sites and 0.84 for the corn sites [64]. In order to compare with modeled flux partitioning, where energy balance is inherently assumed, the measured LE and H were adjusted for closure errors using the approach suggested by Twine et al. [65], which conserves the observed Bowen ratio. For a detailed discussion of the flux tower observations and measurement uncertainty during SMACEX, please refer to Kustas et al. [52] and Prueger et al. [66].

3.2. Hiwater-MUSOEXE Experiment

3.2.1. Description of HiWATER-MUSOEXE Campaign

In order to prove the performance of proposed OSML, we select an additional experimental region with contrasting climatic features and different surface complexities in the Heihe River Basin of northwestern China. Heihe Watershed Allied Telemetry Experimental Research (HiWATER) is designed to be a comprehensive eco-hydrological experiment under the framework of the “Integrated research on the eco-hydrological process of River Basin”. The overall objective of HiWATER is to improve the observability of hydrological and ecological processes, to build a watershed observing system and to enhance the applicability of remote sensing in integrated eco-hydrological studies and water recourse management at the basin scale [67]. HiWATER was a four-year program which formally initialized in May 2012. The Multi-Scale Observation Experiment on Evapotranspiration (MUSOEXE) over Zhangye oasis is part of the HiWATER, and intends to acquire all components and principle parameters related to evapotranspiration processes through multi-scale observations [68]. The underlying surfaces of study area are primarily comprised by oasis cropland, residential areas, sandy Gobi desert (Figure 4). HiWATER-MUSOEXE provides a valuable data source to validate and evaluate the ET models over heterogeneous land surfaces. A detailed description of HiWATER-MUSOEXE can be found in Li et al. [67] and Xu et al. [68].

3.2.2. Site Measurements and Remotely Sensed Data for OSML Input over MUSOEXE

MUSOEXE involved two flux observation matrixes in the middle reach of the Heihe River Basin: one large experimental area and one kernel experimental area (5.5 × 5.5 km²). The large experimental area contained one superstation and four ordinary stations with underlying surfaces of desert (site #19), desert steppe (site #20), Gobi (site #18) and wetland (site #21). The 5.5 × 5.5 km² kernel experimental area is located in the oasis. Overall, 21 eddy covariance system sets, and 21 automatic weather station sets are involved in the experiment. The distribution of flux towers and the land use classifications in the MUSOEXE over Zhangye oasis are shown in Figure 4. Meteorological observations including wind speed, wind direction, air temperature, vapor pressure and atmospheric pressure were sampled at 10-min intervals. Net radiation was measured using a four-component radiometer. Turbulent fluxes were measured using EC system with 10 Hz. The raw data were processed using Edire software and averaged over 30 min. The soil heat fluxes were measured using three heat flux plates located 6 cm below ground surface at each site. The soil heat flux was corrected using the temperature prediction-correction method proposed by Yang and Wang [69] with observed soil moisture profile and soil temperature profile. Sensible and latent heat fluxes were adjusted by forcing the energy balance closure using the measured net solar radiation, soil heat flux, and Bowen ratio [65]. Spatial distributions of meteorological variables, including air temperature, wind speed and relative humidity were interpolated using site observations based on the inverse distance weighting method.

Three scenes from the Advance Spaceborne Thermal Emission and Reflection Radiometer (ASTER) on board Terra acquired on 10 July (DOY 192), 2 August (DOY 215) and 18 August (DOY 231) of 2012 were used in this study. The visible and near-infrared bands of ASTER images have a spatial resolution of 15 m, and thermal infrared bands have a spatial resolution of 90 m. The LST data were retrieved by Li et al. [70] from the thermal infrared bands of ASTER using temperature/emissivity separation (TES) algorithm proposed by Gillespie et al. [71], combined with the water vapor scaling atmospheric correction method [72]. For details concerning the atmospheric correction and the retrieval of LST, please refer to Li et al. [70]. Land surface albedo in the Heihe River Basin in 2012 was retrieved from charge-coupled device (CCD) camera on HJ-1 satellite by Sun et al. [73] in this study. The LST and the land surface albedo are provided by “Heihe Plan Science Data Center, National Natural Science Foundation of China” (http://www.heihedata.org). The albedo data have a spatial resolution of 30 m. For consistency with the thermal infrared band, the visible, near-infrared (NIR) bands of ASTER and albedo data set from HJ-1 were re-sampled to 90-m resolution. Empirical relationships between ASTER-derived vegetation index and LAI, crop height and fractional cover, over the watershed were given by Yang et al. [74]. In this study, these empirical relationships were used to derive regional vegetation parameters for the HiWATER-MUSOEXE region.

4. Results

4.1. Validation at SMACEX Site

The OSML and SEBS models were applied to the SMACEX region using input from satellite-derived and ground-based meteorological observations. The estimation from both OSML and SEBS were firstly validated against tower-based observations, and then the spatial patterns of LE and H generated from both models were inter-compared, which is followed by pixel-based comparison of radiometric-convective resistance estimation from both models.

4.1.1. Validation with Tower-Based Observations

To validate against the observations from the SMACEX tower network, the estimated fluxes were averaged over the upwind source-area (1–2 pixels) for each flux tower [53,54,75]. In addition, for fair comparison, the OSML and SEBS models were evaluated under the same amount of available energy, i.e., driven by the same amount of R_n and G. Validation of energy balance components of R_n and G are shown in Figure 5a,b, while the OSML- and SEBS-derived LE and H estimations are plotted against corresponding EC measurements in Figure 5c,d respectively. The validation statistics are summarized in

Table 2. Summary of validation statistics for validating estimated surface fluxes (W/m²).

Flux	OSML			SEBS
Component	Bias	MAD	RMSD	Bias	MAD	RMSD
Rn	7.3	18.5	24.0	7.3	18.5	24.0
G	−6.9	20.1	23.3	−6.9	20.1	23.3
LE	23.5	39.6	46.5	51.2	74.2	83.2
H	−18.4	30.2	34.5	−41.6	58.0	69.2

.

Results show that the estimated R_n and G agree well with corresponding ground observations. Model bias in the estimated R_n is lower than 10 W/m² and estimated G has a root mean square deviation (RMSD) of 23.3 W/m² and mean absolute deviation (MAD) of 20.1 W/m². The RMSD for H and LE from OSML are 34.5 W/m² and 46.5 W/m², comparatively lower than those from SEBS (RMSD are 69.2 W/m² and 83.2 W/m², respectively). Both models seem to slightly underestimate H and slightly overestimate LE.

4.1.2. Spatial Distribution of Surface Heat Fluxes

The spatial patterns of LE and H estimated by OSML and SEBS based on the Landsat-derived retrievals are shown in Figure 6 and Figure 7 respectively. The spatial patterns derived from OSML and SEBS are quite similar, and the contrast between soybean and cornfields is clear in both OSML- and SEBS-derived H and LE maps. The LE of the cornfields is higher than that in the soybean fields by about 100 W/m². The statistics for regional-scale LE and H estimations are listed in

Table 3. Statistics for regional-scale land surface heat fluxes estimation (W/m²).

			DOY174		DOY182
			Mean	Standard Deviation	Mean	Standard Deviation
LE	OSML	Regional	390.6	83.0	386.8	92.3
		Corn	425.4	59.1	450.5	62.2
		Soybean	313.8	61.0	303.3	56.6
	SEBS	Regional	387.9	107.8	413.4	114.0
		Corn	411.0	104.7	470.2	103.4
		Soybean	324.8	68.3	344.1	76.8
H	OSML	Regional	93.9	55.2	117.6	57.5
		Corn	78.1	39.2	84.4	38.6
		Soybean	146.7	42.9	175.7	38.8
	SEBS	Regional	95.7	83.6	101.9	75.0
		Corn	91.6	87.6	74.9	75.3
		Soybean	135.3	47.8	139.0	53.2

. The mean LE and H derived from the two models are very close at both regional and field scale (i.e., in cornfields and soybean fields).

A pixel-based intercomparison of LE and H estimated by OSML and SEBS was further conducted. Because the employed empirical relationships between crop height and LAI are derived from in-situ observations in cornfields and soybean fields from Anderson et al. [60], this empirical relationship may not be applicable to other land use types. Thus, intercomparison was restricted to the pixels located in the cornfields and soybean fields. The scatter plots of the LE and H derived from OSML against those derived from SEBS in the cornfields and soybean fields are shown in Figure 8. In order to avoid overlapping effect in the scatter plot, we use color shading to indicate pixel density. Generally, OSML-derived LE shows good agreement with SEBS-derived LE. The majority of pixels are closer to 1:1 line, indicating that LE values derived from both models are very similar. Discrepancies in modeled LE are more obvious in pixels with lower LE. The corresponding H values derived from the two models are also presented in Figure 8. The discrepancies are clearly discerned in large H values, indicating that H estimations derived from the two models are discernible in partially-vegetated areas. The performances of different models were highlighted in partially-vegetated landscapes in existing literature. For instance, a similar study by Choi et al. [53] reported significant discrepancies between TSEB- and METEIC-based H estimations in areas with partial vegetation cover (LAI less than 2) in the SMACEX site. In this study, both SEBS and OSML are one-source models, and the difference in the estimated radiometric-convective resistance may be responsible for the discrepancies of H estimation from these two models.

4.1.3. Intercomparison of Radiometric-Convective Resistance Derived from OSML and SEBS

The determination of the radiometric-convective resistance has great impact on the performance of one-source model. In this study, radiometric-convective resistances derived from both SEBS and OSML (denoted as r_{ae_SEBS} and r_{ae_OSML} respectively) are inter-compared to further investigate the difference in performance of these two models. In SEBS, a variable kB⁻¹ approach is applied to estimate radiometric-convective resistance. In OSML, the radiometric-convective resistance was estimated in combination with the WDI derived from the conceptual VFC/LST trapezoid space and the electrical analog formula of H (please see Section 2.2 and Appendix A for details). The statistics for r_{ae_OSML} and r_{ae_SEBS} for cornfields and soybean fields are listed in

Table 4. Statistics for regional-scale r_{ae_OSML} and r_{ae_SEBS} for cornfields and soybean fields (s/m).

		DOY174	DOY182
		Mean	Mean
rae_OSML	Corn	51.8	43.6
	Soybean	52.3	48.5
rae_SEBS	Corn	47.5	56.0
	Soybean	60.4	65.0

. The mean values of r_{ae_OSML} and r_{ae_SEBS} are very similar in both cornfields and soybean fields, which are in alignment with results in Section 4.1.2, that the regional mean of LE and H derived from the two models are very close at both cornfields and soybean fields. To further compare the differences in estimated radiometric-convective resistance from both OSML and SEBS, a pixel-based comparison of r_{ae_OSML} and r_{ae_SEBS} in cornfields and soybean fields on DOY 182 was conducted, as shown in Figure 9. Although different strategies are adopted in SEBS and OSML, r_{ae_OSML} and r_{ae_SEBS} are comparable, which can be observed from the scatter plots of Figure 9. However, Figure 9 also illustrates the differences between OSML- and SEBS-derived r_ae, and this is consistent with the discrepancy in H over partially-vegetated areas observed in Figure 8b,d.

In order to evaluate radiometric-convective resistance derived from SEBS and OSML, the radiometric-convective resistance is calculated from in-situ observations of H (denoted as r_{ae_obs}), and used as an important benchmark here. The evaluation of r_{ae_OSML} and r_{ae_SEBS} for cornfields and soybean fields against observed r_ae is shown in Figure 10. Results show that both SEBS and OSML can provide relatively accurate estimate of radiometric-convective resistance. The RMSD and Bias for r_{ae_SEBS} are 39.0 s/m and 27.9 s/m, whereas the RMSD and Bias for r_{ae_OSML} are noticeably lower (18.8 s/m and 9.6 s/m, respectively). This validation of r_ae lends further credibility to the proposed OSML in surface fluxes estimation.

4.2. Validation at HiWATER-MUSOEXE Site

4.2.1. Validation with Tower-based Observations

Similarly, both OSML and SEBS were applied to the MUSOEXE site using input from satellite retrievals and regional meteorological observations. Consequently, the performance of OSML was evaluated using flux measurements from the observation matrix. Flux measurements were linearly interpolated to match the time of satellite overpass. Figure 11 shows the validation scatter plot for each energy balance component (R_n, G, LE, and H) from OSML and SEBS estimation at the time of satellite overpass. Statistics comparing OSML and SEBS output and flux measurements are summarized in

Table 5. Statistics of validating OSML and SEBS using tower-based observations (W/m²).

Flux	OSML			SEBS
Component	Bias	MAD	RMSD	Bias	MAD	RMSD
Rn	−0.2	20.8	27.8	−1.6	20.3	27.8
G	−23.5	33.3	38.6	−4.5	26.2	42.5
LE	11.8	50.4	67.0	24.5	51.8	67.8
H	13.0	36.5	50.1	−23.3	49.4	57.4

. Results show that the R_n estimations from both models are in good agreement with tower-observed observations, with the absolute mean bias in estimated R_n less than 10 W/m². At HiWATER-MUSOEXE site, the empirical model proposed by Allen et al. [59] and Bastiaanssen et al. [76] was adopted to estimate G, which produced less credible G estimation relative to the estimated R_n, as the scatter plot between model derived G and observations shows a more scattered pattern. Estimated G has a RMSD of 38.6 W/m² and MAD of 33.3 W/m² for OSML. In terms of surface heat flux, SEBS tends to overestimate H, with the MAD and RMSD of H being 49.4 W/m² and 57.4 W/m². Overall, the performance of OSML is comparable to SEBS and to other published studies [74,77,78].

4.2.2. Spatial Patterns of Estimated Regional Land Surface Heat Fluxes

The spatial distributions of H and LE over MUSOEXE site based on OSML and SEBS at satellite overpass time on DOY 192, DOY 215 and DOY 231 are shown in Figure 12 and Figure 13, respectively. The spatial patterns produced by OSML and SEBS are very similar. The contrast between sandy Gobi desert and Zhangye oasis are clear in OSML and SEBS-derived LE and H maps. The Zhangye Oasis, which mainly comprises irrigated farmland, exhibits LE values lager than 350 W/m². On the contrary, the sandy Gobi desert, where precipitation is the only source for soil water, has LE values less than 300 W/m² and H values larger than 100 W/m². The spatial pattern of estimated H and LE on DOY 192 and DOY 215 are quite similar. However, some spatial changes were found on the estimated H and LE map on DOY 231. Specifically, LE increased in the surrounding desert steppe and Gobi desert due to the 2–4 mm of rainfall on DOY 231.

In order to further compare the differences in estimated surface heat fluxes from both OSML and SEBS, a pixel-based comparison of H and LE in the kernel experimental area (the yellow rectangular in the Figure 4) of MUSOEXE was conducted and shown in Figure 14. The fact that most points are concentrated along the 1:1 line indicates that heat fluxes derived from both models are comparable. In addition, the difference between OSML- and SEBS-derived H and LE can also be discerned. Specifically, LE derived from SEBS seems higher than OSML estimations, while H derived from SEBS is lower than OSML. The statistics comparing OSML and SEBS are summarized in

Table 6. Statistics for comparing OSML- and SEBS-derived H and LE for the kernel experimental area in MUSOEXE (W/m²).

	DOY 192		DOY 215		DOY 231
	LE	H	LE	H	LE	H
Bias	24.4	−0.5	36.0	−13.1	28.1	−21.1
MAD	33.8	12.9	47.0	16.8	43.4	21.4
RMSD	38.3	21.0	49.5	19.6	47.1	22.6

, which further indicate H and LE derived from two models are similar.

5. Discussion

5.1. Advantages of OSML over Typical One-Source Models

Due to the discrepancy between LST and T_aero, the determination of radiometric-convective resistance poses a major challenge for regional heat flux estimation based on remotely sensed LST. In this study, a One-Source Model for Land (OSML) is proposed to address this issue. The performance of OSML is comparable to the typical one-source model of SEBS when validated with the SMACEX data and HiWATER-MUSOEXE data, both of which have been widely used to validate and evaluate the ET models [4,34,53,54,64,79,80,81,82]. The prominent merits of OSML are as follows: (1) In OSML, the radiometric-convective resistance is estimated by combing the WDI derived from the conceptual VFC/LST trapezoid space and the electrical analog formula of H. This manner of physically deriving radiometric-convective resistance eliminates the need for an empirical extra resistance, which has long been a difficulty in land surface heat fluxes estimation with remotely sensed LST. (2) Information on wind velocity and aerodynamic surface characteristics is not required in OSML. In most one-source and two-source models such as SEBS and TSEB, wind speed is an indispensable model input. However, variation of wind velocity at regional scale is difficult to reproduce using in-situ observations [17,32,33]. In addition, the aerodynamic surface characteristics are often determined empirically based on land use types or crop height, which will inevitably cause considerable uncertainties. The independence from information of wind speed and aerodynamic surface characteristics has made OSML more suitable for application at regional scale with heterogeneous surface.

5.2. Discussion on Different Approaches to Mediate Difference between LST and T_aero

As stated in Section 1, to account for the difference between LST and T_aero, several approaches have been proposed to adjust aerodynamic resistance so that H can be estimated using LST directly. For instance, Chehbouni et al. [20] proposed a β function (β_ch) and Matsushima [22] proposed a α function to interpret the relationship between LST and T_aero. The related adjustment parameters are often estimated based on its empirical relationship with LAI. Most of the works were based on the adjustments of the kB⁻¹ parameter. Among them, a variable kB⁻¹ approach is adopted in SEBS model. However, this approach is found to be unstable under sparse steppe [28] and sparsely vegetated area [77]. The work by Lhomme and Montey [21] and Zhao et al. [83] indicate that a universal parameterization of kB⁻¹ seems impossible.

As an alternative to the above-mentioned approaches, radiometric-convective resistance is used in lieu of aerodynamic resistance for H estimation with radiometric surface temperature. In this study, the evaluation of r_{ae_SEBS} and r_{ae_OSML} against r_{ae_obs} suggests that OSML provides a novel and reliable way to parameterize radiometric-convective resistance (Figure 10). In comparison, the variable kB⁻¹ approach used to determine radiometric-convective resistance in SEBS requires the parameterization of crop height, surface roughness and displacement height. In addition, wind speed, which is difficult to obtain at regional scale, is an indispensable variable for SEBS model. In contrast to SEBS model, OSML does not need to quantify aerodynamic related parameters, which will greatly promote the application of surface heat fluxes estimation at regional scale.

5.3. Operationality of OSML

Satellite remote sensing is a promising technique for obtaining regional land surface heat fluxes and this information is important for precise quantification of the water balance for water resources planning and management, optimizing crop production, identifying crop stress and drought impacts [3,8,9,58,59]. In this study, the regional application of OSML employs Landsat images, whereas the Moderate Resolution Imaging Spectroradiometer (MODIS, 1 km resolution) or geostationary meteorological satellite images from Advanced Very High Resolution Radiometer (AVHRR, 5 km resolution) would also provide sufficient data source for routinely monitoring surface heat fluxes using OSML. Besides, regional meteorological variables required for routinely application of OSML are also readily available from China Meteorological Administration. The real-time or near-real-time reanalysis products from a data assimilation system, such as NASA Global Modeling and Assimilation Office’s (GMAO), Modern Era Retrospective Analysis for Research and Applications (MERRA) data set [84], the Global Land Data Assimilation System data set [85], and China Meteorological Administration (CMA) Land Data Assimilation System (CLDAS) data set [86], would make it possible to routinely apply OSML at regional scale for ET monitor/prediction. For longer time applications, developing a reliable method for extrapolating daily values to longer time scales becomes imperative.

6. Conclusions

In this study, a One-Source Model for Land (OSML) is proposed for estimation of regional-scale surface fluxes. In OSML, the conceptual VFC/LST trapezoid model and the electrical analog formula of H are combined to analytically resolve the radiometric-convective resistance, and thus circumvent the difficulty of parameterizing the empirical “extra resistance” in conventional one-source energy balance models.

To evaluate the performance of the proposed method, the OSML was applied over two regions with contrasting climatic features, namely, the SMACEX and HiWATER-MUSOEXE. Compared with tower observations, the RMSD of H and LE from OSML are 34.5 W/m² and 46.5 W/m² at SMACEX site and 50.1 W/m² and 67.0 W/m² at HiWATER-MUSOEXE. Overall, the performance of OSML is very comparable to typical one-source models such as SEBS and other published studies. A further intercomparison of radiometric-convective resistance derived from both OSML and SEBS indicates that the proposed novel method can produce reliable radiometric-convective resistance and thus accurately estimate surface fluxes, without the requirement of information on wind speed and specification of surface characteristics. Overall, the study suggests that the OSML is promising and can be used to estimate regional surface heat fluxes over heterogeneous surface.