The Two-layer Surface Energy Balance Parameterization Scheme ( TSEBPS ) for estimation of land surface heat fluxes

A Two-layer Surface Energy Balance Parameterization Scheme (TSEBPS) is proposed for the estimation of surface heat fluxes using Thermal Infrared (TIR) data over sparsely vegetated surfaces. TSEBPS is based on the theory of the classical two-layer energy balance model, as well as a set of new formulations derived from assumption of the energy balance at limiting cases. Two experimental data sets are used to assess the reliabilities of TSEBPS. Based on these case studies, TSEBPS has proven to be capable of estimating heat fluxes at vegetation surfaces with acceptable accuracy. The uncertainties in the estimated heat fluxes are comparable to in-situ measurement uncertainties.


Introduction
Land surface Evapotranspiration (ET) is one of the most important components in the water cycle between the earth and atmosphere, and plays a very important role in the atmosphere, hydrosphere, and biosphere of the planet.It is an urgent task to understand the evapotranspiration process over different surface types and conditions in agriculture, hydrogeology, forest, and ecology for the purpose of using water resources properly.Additionally, land surface evapotranspiration is a key parameter in the synoptic and climatic phenomenon because of the heat and moment transfer processes in association with evapotranspiration.Studies (Dickinson, 1984;Avissar, 1998) on climate models and general circulation models (GCMs) have found that the climate is sensitive to the change of land surface evapotranspiration.At present, remote sensing may be the only efficient technical way that Correspondence to: X. Xin (xin xzh@sohu.com)can be used to monitor surface evapotranspiration on the regional scale (Mu et al., 2007;Stisen et al., 2008).Spatial and temporal distributions of the key state variables of the land surface energy balance can be provided by remote sensing, and can be used to estimate surface evapotranspiration.The data of mid-low resolution meteorology and the land resource satellite can cover large areas of the land surface and can observe repeatedly in short periods, which is useful for the research in the drought monitoring, climate changes, water resource management, and so on.
Generally, surface evapotranspiration (i.e.latent heat flux LE) is estimated as the residual term of surface energy balance equation.Remotely sensed data have been used successfully over the past years to estimate the surface net radiation and the soil heat flux (hence available energy) from combined visible, near infrared and thermal infrared data (Norman et al., 1995;Liang et al., 2000;Jacobs et al., 2000;Ma et al., 2002;Ma, 2003).Therefore, the primary focus has been the determination of the sensible heat flux based on the spatially distributed surface temperature fields.The turbulent heat fluxes models to estimate the sensible heat flux can be categorized into two groups, single-source models and dualsource models, according to whether or not the model separates the foliage and the substrate soil.In the single-source models, a so called "excess" resistance or parameter kB −1 is used to account for the difference between the remotely sensed radiative surface temperature T r and the aerodynamic temperature T 0 (Moran et al., 1989;Kustas, 1990).The difference between T 0 and T r depends on a number of factors within the Soil-Plant-Atmosphere Continuum (SPAC) as well as the viewing condition of the Thermal Infrared (TIR) sensor.Therefore, it is very difficult to find out a robust relationship that takes all these factors into account (Choudhury et al., 1986;Troufleau et al., 1997;Chehbouni et al., 2001).Many authors (Blyth et al., 1995;Verhoef et al., 1997;  Troufleau et al., 1997;Kustas et al., 1999;Massman, 1999) have examined the features of the kB −1 parameter.This parameter is a complex function of canopy structure, water stress and environment factors, and it is too variable to provide a universal solution for estimating the sensible heat flux using single-angle radiative surface temperature.This problem can be circumvented to some extent by using the dualsource models.In this type of models, the heat fluxes of the components (foliage and soil) are simulated individually, and the aerodynamic temperature is analytically expressed in terms of the component temperatures and a set of resistances, as described in the two-layer model proposed by Shuttleworth and Wallace (1985) and revised by Shuttleworth and Gurney (1990).This is very important for sparsely vegetated surfaces, because in this circumstance the contribution of soil surface cannot be neglected.Otherwise, the bias of the estimated surface heat fluxes can be significant.
Even though the advantage of the dual-source models in physics has been recognized by the scientific community, the most widely used methods in applications are still based on the assumption of the single source of the surface heat fluxes.This results from such a fact that the use of the twolayer model for operational purpose requires component surface temperatures (i.e.soil and vegetation), which is still not available from regular observations and retrieval of the most space-borne remote sensors.Studies of applying the twolayer model with traditional single-angle TIR data have been reported since the model was proposed (Norman et al., 1995;Jupp et al., 1998).Usually, this is achieved by simplification of the model or adding an empirical relationship in the model, which decreases the modeling accuracy or limits universal application.
In this study, we have developed a physics-based Two-layer Surface Energy Balance Parameterization Scheme (TSEBPS) for estimation of land surface heat fluxes.We combined the two-layer model developed by Shuttleworth and Wallace (1985) with techniques of handling limiting cases as shown in Su (2002) and Norman et al. (1995) to derive the Component Temperature Difference (CTD) under several extreme soil moisture states.
Additionally, a directional thermal radiative transfer model is used to simulate the radiative surface temperature at these states.Then an index is developed using the observed surface temperature and the simulated temperature at the extreme states.This index is then used to calculate the actual sensible and latent heat fluxes of the foliage and soil surface.

TSEBPS (Two-layer Surface Energy Balance
Parameterization Scheme)

The Two-layer Surface Energy Balance model
The classical two-layer model by Shuttleworth and Wallace (1985) founded the theory basis for this study (Fig. 1).
The surface energy balance is commonly written as Where R n is the net radiation, G is the soil heat flux, H is the sensible heat flux, and LE is the latent heat flux (L is the latent heat of vaporization and E is the actual evapotranspiration).The net radiation of the surface (R n ) can be calculated from the equation: Where S d are solar irradiation, α surface albedo, ε s surface emissivity, L d downward atmosphere long wave radiation, and L u surface emitted long wave radiation.G Can be calculated with method used by Su (2002): Where, s =0.315 and c =0.05, and f c fractional canopy coverage.
The budget of the net radiation between soil and the canopy can be calculated using the Beer's law: Where R ns and R nv are the net radiation of soil and the canopy, and b(θ) is the gap frequency of the canopy written as Where, θ is the solar zenith angle, LAI leaf area index of the canopy, and G(θ) projection coefficient of the leaves which is related to the Leaf Angle Distribution (LAD).The energy balance of the soil is written as: The energy balance of the canopy is written as: Hydrol.Earth Syst.Sci., 14, 491-504, 2010 www.hydrol-earth-syst-sci.net/14/491/2010/ The basic principle underlying two-layer models is that the two sources of water vapor and heat are superimposed and hence heat and water vapor enter or leave the bottom layer only via the top one.The total flux of sensible heat emanating from the whole surface is the sum of the fluxes emanating from each layer (here soil and vegetation).So there is where, ρ is the air density (kg m −3 ), C p the specific heat of air at constant pressure (J kg −1 K −1 ), T 0 the aerodynamic temperature (K) defined as the extrapolation of the air temperature profile down to the apparent source/sink of heat within the canopy, T a air temperature (K) at the reference height, and r aa the aerodynamic resistance (s m −1 ) for heat transfer.H s and H v are soil and vegetation sensible heat fluxes, respectively, which can be expressed according to the gradient-diffusion hypothesis as Where, T s and T v are soil and vegetation temperature, respectively, r a s the aerodynamic resistance between soil and the source height in the canopy, and r a v the bulk boundarylayer resistance of the vegetation.The transfer of the latent heat flux in the canopy can also be expressed similarly as: where, γ is the psychometric constant (kPa K −1 ), e 0 the aerodynamic vapor pressure of the surface, e a vapor of the atmosphere, LE s and LE v soil and vegetation latent heat fluxes respectively, e(T s ) and e * (T v ) vapor pressure of soil surface and the saturation vapor pressure in leaf stomata respectively, r ss , and r st soil surface resistance and leaf stomata resistance respectively.Aerodynamic resistance r aa is formulated using the stability correction method by Choudhury (1989): Where r a0 is the aerodynamic resistance in the neutral atmosphere condition: Where u is the wind speed at the reference height z, and k von Karman's constant.The corrective term φ is calculated with: Where g is acceleration due to gravity (ms −2 ).The zero plane displacement height d and the roughness length for momentum z 0 can be determined following Choudhury and Monteith (1988), who fitted simple functions to the curves obtained by Shaw and Pereira (1982) from the second-order closure theory: (16) Where, c d is the mean drag coefficient assumed to be uniform within the canopy (0.2), and z 0s the roughness length of the substrate.For bare soil, z 0s is taken as 0.01 m.The formulations for resistances r a s and r a v proposed by Choudhury and Monteith (1988) and Shuttleworth and Gurney (1990) are used here: Where w is the leaf width, u(h) the wind speed at canopy height h, α 0 and α w two constant coefficients equal to 0.005 (ms −1/2 ) and 2.5 (dimensionless), respectively.The value of eddy diffusivity at canopy height K (h) is determined with K (h)=ku * (h−d).

Parameterization scheme based on limiting cases
Figure 2 gives the flow chart of the parameterization.First of all, the limiting cases of soil moisture in the Soil-Plant-Atmosphere Continuum (SPAC) are defined, which are drylimit, wet-limit, and transition-state.The definitions of the dry-and wet-limit are similar to those in SEBS (Su, 2002), but differ in processing soil and foliage components individually.The transition-state occurs when the surface soil layer is dry and the root zone soil is still wet, which is understandable and predictable in natural vegetation because the drying-off process after a rainfall or irrigation event starts from the surface.Then the component temperature difference (CTD, i.e., T s −T v ) at the limiting cases is derived based on the following assumptions.
Under the dry-limit, the latent heat (or the evaporation and transpiration) becomes zero due to the limitation of soil  moisture and the sensible heat flux is at its maximum value.From Eqs. (1), ( 7) and (8), it follows, and The CTD under this case can be derived from Eq. ( 10).
The aerodynamic surface temperature at dry-limit T 0,dry can also be calculated from Eq. ( 9) based on above assumption.Hence, the soil and foliage temperatures under this case T s,dry and T v,dry can be calculated using δT dry and T 0,dry .Under the wet-limit, where the evaporation and transpiration take place at potential rates (i.e. the evaporation and transpiration is limited only by the energy available under the given surface and atmospheric conditions), the sensible heat flux takes its minimum value.
The aerodynamic surface temperature at wet-limit T 0,wet can be calculated from Eq. ( 9) based on above assumption.The component temperature difference between soil and foliage can be derived based on the P-M type equation of soil and the canopy and assuming the soil surface resistance and the stomata resistance are zero, we have where δT wet is CTD under the wet-limit, wet is the slope of the saturation vapor pressure versus the temperature, and γ is psychrometric constant.Hence, the soil and foliage temperatures under this case T s,wet and T v,wet can be calculated using δT wet and T 0,wet .
Under the transition-state, where the evaporation becomes zero due to the limitation of surface soil moisture, and the transpiration is limited only by the energy available (i.e., root zone soil moisture is still at wet-limiting).So there is: and the transpiration is simulated using Priestly-Taylor equation.
where Priestly-Taylor constant a=2.0 according to Kustas et al. (1999), f g is fraction of green leaves in the canopy.So the aerodynamic surface temperature T 0,trans and foliage temperature T v,trans under this case can be calculated using Eqs.( 9) and ( 10), and the soil temperature T s,trans under this case is derived using T 0,trans and T v,trans .Based on the above assumptions and calculations, we have the aerodynamic surface temperature under the limiting cases, T 0,dry , T 0,wet , and T 0,trans , and the soil and foliage temperatures under the limiting cases, T s,dry , T v,dry , T s,trans , T v,trans , and T s,wet , T v,wet .So we also have the sensible and latent heat fluxes of the soil and foliage under the limiting cases, H s,dry , H v,dry , LE s,dry , LE v,dry , H s,trans , H v,trans , LE s,trans , LE v,trans , and H s,wet , H v,wet , LE s,wet , LE v,wet based on Eq. ( 10).
The next step is to derive the actual sensible and latent heat fluxes of the soil and foliage using an interpolation method from the limiting cases.We assume that the dryand wet-limit cases set reasonable boundaries of the surface heat balance under limiting conditions, and the transitionstate gives a key spot where dramatic changes of the budget of sensible and latent heat of the canopy take place (i.e., transpiration is at its maximum value and evaporation decreases between wet-limit and transition-state, and evaporation is zero and transpiration decreases between transitionstate and dry-limit).Increasing or decreasing the soil and foliage heat fluxes can bring about changes in the temperatures of the soil and foliage, which can result in canopy surface temperature changes.We have derived the component temperatures under the limiting-cases, from which we simulated the radiometric surface temperature under the limiting cases, T r,dry , T r,wet , andT r,trans using a directional thermal infrared radiative transfer model of the canopy.In this study, the model proposed by Franc ¸ois (1997) was used to simulate directional radiometric surface temperatures.In the simulation, the observing zenith angle takes the actual angle in the field measurement of T r , and the soil and foliage emissivity takes the value of 0.94 and 0.98 following Franc ¸ois (1997) andFranc ¸ois (2002).So the actual heat fluxes can be derived based on the comparison between the actual surface temperature and the simulated surface temperature under the limiting-cases.
Comparison between the measured radiometric surface temperature and the simulated surface temperature under the limiting cases can give a clue of the status of soil moisture, i.e., higher temperature than that under the transition state hints limitation of soil moisture on evaporation, and lower temperature than that under the transition state may indicate relatively better soil moisture condition in the canopy.The derivation of the actual heat fluxes is: (1) If T r,wet <T r <T r,trans , transpiration is at its maximum value and evaporation decreases with increasing surface temperature, we have: where x is an index build from radiometric surface temperatures: The sensible heat flux of soil and foliage is then derived as the residual of the energy balance equation of the soil and foliage.
(2) If T r,trans <T r <T r,dry , soil sensible heat flux is at its maximum value (evaporation is zero) and foliage sensible heat flux increases with increasing surface temperature, we have: where y is an index build from radiometric surface temperatures: The latent heat flux of soil and foliage is then derived as the residual of the energy balance equation of the soil and foliage.
The indices x and y are used to measure the relative distance of the actual radiometric surface temperatures from the virtual radiometric surface temperatures under the limiting cases.The coefficient n is used to account for the non-linear effect of the heat fluxes changing with the relative change of the surface temperature.Here we take the value of n=0.25 and it shows that the result is not sensitive to this coefficient.
(3) If an unexpected situation happens, such as T r >T r,dry or T r <T r,wet , which may result from the errors of the measurements, simulations and assumptions, the heat fluxes under the limiting cases are used for the actual heat fluxes.

Data
Two sets of in-situ data were used for evaluation of TSEBPS: (1) the data set from the "Quantitative Remote Sensing theory and application for Land Surface Parameters (QRSLSP)" project at Shunyi, Beijing, China, 2001, and (2) the data set from the "Watershed Allied Telemetry Experiment Research (WATER)" project in the Heihe River Basin, Gansu, China, 2008.

Winter wheat in Beijing
The winter wheat dataset was obtained during the "Quantitative Remote Sensing theory and application for Land Surface Parameters (QRSLSP)" campaign that was carried out in North China in April 2001.The main concern of this experiment was for quantitative remote sensing applications in agriculture.The winter wheat fields located in Shunyi district, north of Beijing ( 116• 34 E, 40 • 12 N) were selected as the chief observation target.The winter wheat with row structure and regular irrigation is one of the main agricultural crops in North China, and usually the growing period after the winter starts from the end of March through the beginning of April.The experiment was carried out in April in order to obtain the in-situ data during the rapid growing period of the winter wheat.There are three observation sites, NW3, NW4 and NW5 that are adjacent from south to north, with different planting and management measures, such as wheat cultivar, sowing date, irrigation/fertilization date and amount due to the fields belonging to different farmers, which resulted in different surface conditions among the three sites especially the soil moisture.During the experiment period, soil moisture condition was the best in NW4 and the worst in NW5, which resulted in evident difference in heat fluxes and surface temperature between the fields.Turbulent heat fluxes and meteorological data were measured with Bowen-Ratio (BR) system and Automatic Weather Station (AWS) at the 3 sites, respectively (see Table 1).The interchange of high-and low-layer measurements takes place for every 10-min for sites NW3 and NW4, and 5min for site NW5, from which 20-min (NW3 and NW4)/10min (NW5) average turbulent fluxes (H and LE) were computed in order to eliminate the discrepancy of equipments at the two sides of the system.10-min averages of net radiation and soil heat flux were stored.The measured soil heat flux is the value at the 5 cm under the surface for the all sites in this study, and was corrected to the surface by the method of integration using the gradient of soil temperature and the soil heat flux (Liebethal et al., 2005).In addition, 10-min averaged ancillary meteorological data, such as air temperature, relative humidity, and wind speed were also recorded.10min average surface brightness temperature was measured and recorded by TIR radiometers, from which the radiative surface temperature was obtained by correction of atmospheric effect and emissivity (Olioso et al., 1996).Hence, every 20-min (NW3 and NW4)/10-min (NW5) averaged heat fluxes, net radiation, soil heat flux, meteorological data, and surface temperature during daytime (when both sensible and latent heat fluxes are positive) were collected as a group of data, and regarded as a sample (see Table 2).The period of available data of the 3 sites are different due to the different beginning/ending time of TIR observation.
As a necessary input for the model, canopy structure data (including Leaf Area Index -LAI, canopy height, leaf shape, and row width and space) were also measured manually by a specific team at the 3 sites regularly during the experiment.
So the winter wheat dataset contains 3 sub-datasets, which represent different soil moisture condition as well as different vegetation density as shown in Table 2.The 3 sub-datasets are used independently to evaluate TSEBPS.More detailed information about the experiment can be found in Liu et al. (2002) for the interested.

Maize in Gansu
The maize dataset was obtained during "Watershed Allied Telemetry Experimental Research (WATER)" project carried out in Heihe River Basin of Gansu province, Northwest China from May to July 2008 (Li et al., 2009).The main concern of this experiment was to study hydrology and ecology processes using remote sensing techniques, therefore evapotranspiration is one of the main concerns in this project.Heihe River Basin of Gansu province is located in the arid/semi-arid region in the northwest of China, where the agricultural and natural ecosystems suffer from deficit of precipitation frequently.The agriculture is supported mostly by Hydrol.Earth Syst.Sci., 14, 491-504, 2010 www.hydrol-earth-syst-sci.net/14/491/2010/ the irrigation system, which takes the melted snow/ice water from the upper-stream Qilian mountain area to the flat middle-and lower-stream oasis.
The site Yingke (YK) is located in the artificial oasis to the south of Zhangye city (100 • 24 E, 38 • 51 N), where the main crop is maize with row structure and regular irrigation.The turbulent heat fluxes and meteorological data were measured with Eddy-Covariance system (EC) and Automatic Weather Station (AWS).Half-hourly averaged turbulent fluxes (H and LE) were computed, while 10-min averages of net radiation and soil heat flux were stored.The measured soil heat flux is the value at the 5cm under the surface for the all sites in this study, and was corrected to the surface by the method of integration using the gradient of soil temperature and the soil heat flux (Liebethal et al., 2005).In addition, 10-min average ancillary meteorological data, such as air temperature, relative humidity, and wind speed were also recorded.About 80% energy closure ratio was found in the EC data.Since the two-layer model requires energy conservation, closure in the flux measurements was enforced through a Bowen-ratio method; that is, Bowen-ratio was calculated using H and LE of the EC measurements, and then H BR and LE BR were recalculated with Bowen-ratio method using net radiation and soil heat flux.10-min average surface brightness temperature was measured and recorded by TIR radiometers, from which the radiative surface temperature was obtained by correction of atmospheric effect and emissivity (Olioso et al., 1996).Hence, every 30-min averaged heat fluxes, net radiation, soil heat flux, meteorological data, and surface temperature during daytime (when both sensible and latent heat fluxes are positive) were collected as a group of data, and regarded as a sample (see Table 2).As a necessary input for the surface models, canopy structure data (including leaf area index -LAI, canopy height, leaf shape, and row width and space) were measured manually from 21 May to 15 July throughout the whole growing period before tasseling stage of maize.
Unlike the field campaign of QRSLSP, the experiment of the WATER project had lasted for several months.The data collected during the experiment covers the main growing period of maize, which allows us to evaluate TSEBPS with data of different vegetation coverage states, i.e., from very sparse vegetation at the beginning (LAI<0.5), to very dense vegetation at the end (LAI>5).In order to evaluate the performance of TSEBPS at different canopy coverage, the dataset of maize was separated into 3 subsets according to LAI; that is YK-sparse for the data when LAI<1.0, YKmedium for 1.0<LAI<3.0,and YK-dense for LAI>3.0.
Table 1 gives the brief information about the turbulent fluxes and TIR radiometric measurements.Table 2 lists the datasets or subsets that are used in the evaluation.In summary, the number of data points is mainly decided by (1) the availability of the observation (because of discontinuity of observation), (2) temporal average of data, (3) processing and quality control of BR and EC data, (4) the data number of daytime (because only the data during daytime when both sensible and latent heat fluxes are positive were used here).

Results
The accuracy of TSEBPS will be assessed using the datasets listed in Table 2. Radiative surface temperature as well as ancillary meteorology and canopy structure data were input to the TSEBPS, and the sensible and latent heat fluxes are estimated as discussed previously.All other input variables are measured including net radiation and soil heat flux.The difference between estimation and measurement of the sensible and latent heat fluxes will be analyzed for each of the datasets.

Results of the winter wheat datasets
The canopy sensible and latent heat fluxes predicted versus the measured values of winter wheat sites are shown in Fig. 3. On the whole, TSEBPS estimated heat fluxes agree very well with the field measurements over winter wheat canopies.The performance of TSEBPS at the 3 sites is very close besides the difference in the magnitude of sensible and latent heat fluxes, which can be explained to some extent by the surface condition of the fields.As we have mentioned before, the canopy density and soil moisture condition are different at the 3 sites (Table 2), which resulted in different magnitude of sensible and latent heat fluxes (therefore the Bowen ratio).The average value of available energy (net radiation minus soil heat flux) for site NW3, NW4 and NW5 is 324.5, X. Xin and Q. Liu: TSEBPS for estimation of land surface heat fluxes 331.4 and 205.2 Wm −2 , respectively.The average value of measured sensible heat flux for the 3 sites is 100.5, 55.4 and 73.4 Wm −2 , and the latent heat flux is 224.0, 276.0 and 107.2 Wm −2 , respectively.For the latent heat flux, the best agreement appears at NW4, and followed by NW3 and NW5, and all of the predictions are within acceptable accuracy.The data points are scattered closely to the 1:1 line and the bias is confined mostly to within around 50 Wm −2 , indicating good agreement with measured values.There is no obvious trend of overestimate or underestimate of the heat fluxes.
Table 3 to Table 5 show the error statistics of the predicted heat fluxes.Root-Mean-Squared-Error (RMSE), Mean-Absolute-Difference (MAD) and Mean-Absolute-Percentage-Difference (MAPD) are shown in the tables.RMSE of the 3 sites are all within 35 Wm −2 and MAD within 30 Wm −2 , which means that the predicted heat fluxes agree well with the field heat fluxes observation.Mean and standard deviation of the predicted heat fluxes compare very well with those measured as shown in Table 3 to Table 5.The best agreement is found at NW4 dataset, where both mean and standard deviation of predicted sensible and latent heat fluxes are very close to the measurements.The discrepancy between measurements and simulation is within the uncertainty of turbulent heat fluxes measurements.Coefficients of determination (R 2 ) for sensible and latent heat fluxes are high at the three sites, indicating TSEBPS can predict heat   fluxes with high accuracy.The highest and lowest R 2 of the predicted latent heat flux appear at site NW4 and NW5, respectively.
In order to investigate the bias of TSEBPS-estimated LE, we compared the relationship between the bias and input parameters and found that the surface temperature gradient (surface temperature minus air temperature) is the mostly related factor with the bias as shown in Fig. 4. We can see that the temperature gradient is mostly under 2 K at NW4, and the bias of estimated LE is also small, mostly within ±20 Wm −2 .At point No. 8 (12:00, 13-April), the temperature gradient is the largest (about 8 K), and the bias of estimated LE is also the largest (about −60 Wm −2 ).At NW5, the temperature gradient is much higher than that of NW4 (mostly between 5∼20 K), and the bias of estimated LE is also larger than that of NW4 (mostly within ±50 Wm −2 ).On the whole, the trend of bias is opposite to that of temperature gradient.Similar to NW4, the points with largest bias (LE was much underestimated in Fig. 3) also have very large temperature gradient.
We also investigated the correlation between the bias of TSEBPS-estimated LE and wind speed.It can be seen from Fig. 4 that there is no obvious trend in the correlation for site NW4 and NW5.Generally, wind speed is negatively correlated with the resistances for the transfer of heat in the canopy-atmosphere system, which means that higher wind speed will result in higher sensible heat flux and lower latent heat flux if we employ a simple single-layer model to calculate sensible heat flux and derive latent heat flux using residual method.In this study, however, as we employed an interpolation method to calculate the sensible (or latent) heat flux and derive latent (or sensible) heat flux as the residual of the energy balance equation, the impact of wind speed on the bias of TSEBPS-estimated LE is not that straight forward.From Eq. ( 27) through Eq. ( 30) we can infer that the bias of TSEBPS-estimated LE is much correlated with surface temperature gradient than wind speed.Wind speed can influence component heat fluxes at the limiting cases, but its influences might counteract each other in Eqs. ( 27) and (29).Meanwhile, wind speed can influence the surface temperature gradient (see Fig. 4), which in turn will propagate to the bias of TSEBPS-estimated LE.

Results of the maize dataset
The canopy sensible and latent heat fluxes predicted versus the measured values are shown in Fig. 5. Similar to the winter wheat dataset, the estimated sensible and latent heat fluxes agree very well with the measurement.Table 6 shows the error statistics of the predicted heat fluxes.The average value of available energy (net radiation minus soil heat flux), sensible and latent heat fluxes is 335.4,73.4 and 262.0 Wm −2 , respectively.RMSE and MAPD of the estimated latent heat flux are low and the coefficient of determination (R 2 ) is very high, which means that the TSEBPS-estimated latent heat flux with TIR measurements can reach high accuracy.Mean and standard deviation of the predicted heat fluxes compare very well with those measured as shown in Table 6.
In order to investigate the performance of TSEBPS at different vegetation coverage conditions, the error statistics are recalculated separately for the 3 subsets of the maize according to Table 2.The results are shown in Table 7, from which we can see that there is no evident difference in the R 2 between the subsets, but the RMSE shows much more variability between the subsets, i.e., RMSE increases with increasing LAI.On the other hand, MAPD decreases with increasing LAI.Comparison of mean and standard deviation shows that datasets of medium and dense canopy have larger bias than that of sparse canopy.However, the difference between the subsets is not evident, and the performance of TSEBPS is stable from very sparse to very dense canopies.It means that TSEBPS can estimate heat fluxes accurately above surfaces with different density of vegetation.
The turbulent heat fluxes were measured by Bowen-ratio system in the winter wheat sites and eddy-covariance system in the maize site.Both techniques are popular in experiments.In this study, EC data was processed to meet the energy balance with a Bowen-ratio method (Twine et al., 2000), and BR data was also processed with quality control.Nevertheless, it is hard to compare the different measurement techniques based on the present datasets and give a conclusion about the uncertainties of the measurements in this study.Fortunately, some useful information can be found in the references that analyzed the variation of flux estimation by various micrometeorological techniques based on the datasets obtained in other experiment projects, such as Monsoon'90, FIFE, and ChinaFLUX (Norman et al., 1995;Twine et al., 2000;Massman et al., 2002;Yu et al., 2006).According to the references and other studies that compare model predicted flux with in-situ measurements (e.g., Timmermans et al., 2007), uncertainties of fluxes are about 25∼50 Wm −2 for H and LE measured by EC technique, and about 20% for LE measured by BR technique.The errors of TSEBPSestimated heat fluxes are of similar magnitude with the uncertainties in the measurements, which means that TSEBPS is able to predict surface heat fluxes with acceptable accuracy.

Error analysis
According to the flow chart of TSEBPS (Fig. 2), the actual heat fluxes are derived from the heat fluxes of the limiting cases with an interpolating method.So the error of TSEBPSestimated heat fluxes comes from these two aspects, i.e., the heat fluxes of the limiting cases and the interpolating methods.The sensitivity of the estimated heat flux to the error of the heat flux at the limiting cases is described by the following way.
where Y represents the derived actual heat flux, and Y i the heat flux at the limiting cases (i.e., wet-and dry-limits, and transition state).From Eqs. ( 27) and ( 29), we can see that the non-linear interpolation takes place for soil latent heat flux when T r,wet <T r <T r,trans , and for foliage sensible heat flux when T r,trans <T r <T r,dry .And at other cases, the interpolation is linear.Sensitivity to the error of LE s,wet in Eq. ( 27) and the error of H v,dry in Eq. ( 29) can be expressed in a same way: where A represents LE s,trans /LE s,wet and p for x for Eq. ( 27), and A represents H v,trans /H v,dry and p for y for Eq. ( 29).
According to the assumption of TSEBPS (Eqs.25 and 26), A equals to 0 or is very close to 0 (no negative value of the heat fluxes is allowed in the calculation), which results in that the sensitivity to the error of LE s,wet and H v,dry is nearly ±10%.It means that the error of component heat fluxes at the dry-and wet-limiting cases is propagated to the estimated heat fluxes in a linear way.
Sensitivity to the error of LE s,trans in Eq. ( 27) and the error of H v,trans in Eq. ( 29) can also be expressed in a same way: Because A equals to 0 or is very close to 0, the sensitivity to the error of LE s,trans and H v,trans is very small and can be regarded as 0. It means that the error of component heat fluxes at the transition state has no obvious influence on the estimated heat fluxes.Sensitivity to the error of p (x in Eq. 27 and y in Eq. 29) can be expressed as: Because A equals to 0 or is very close to 0, the sensitivity to the error of p mainly varies with p.The magnitude of p is within the range of [0, 1].When p is close to 0, the sensitivity is small, and when p is close to 1, the sensitivity becomes relatively larger.And the sign of the error in the estimated heat fluxes is opposite to that of p.In our datasets, the average value of p is about 0.5∼0.6,which leads to about ±10∼20% error in the estimated heat fluxes for ±10% of error in p.
From above analysis we can see that ±10% error in the component heat fluxes at the wet-and dry-limiting cases will result in about ±10% error in TSEBPS-estimated heat fluxes, and the error in the component heat fluxes at the transition state will result in no obvious error in TSEBPS-estimated heat fluxes.The component heat fluxes at the limiting cases are calculated using Eq. ( 10) with the aerodynamic temperature and component temperatures, which are calculated based on the assumptions of the limiting cases.In this study, the assumptions and calculations are physicsbased and the error in the estimated component heat fluxes is regarded within acceptable range.
On the other hand, the error in the simulated surface temperature at the limiting cases has obvious influence on the results.The error of p comes from the error of TIR observation, as well as the error of the simulated surface temperature at the limiting cases.In our study, a directional canopy TIR radiation transfer model by Franc ¸ois (1997) is used to simulate the surface temperature at the limiting cases.This model is of reasonable physics-basis and has performed well in the experimental study in the reference.In their study, the error of the simulated temperature is relatively small and acceptable.In this study, we believe that the simulated temperature is of good quality and comparable to the field TIR observation.Furthermore, from Eqs. ( 28) and (30) we can see that the error in x and y can be relatively small because the index is constructed by the difference between the temperatures, which means that the error of the temperature can wipe one another out.
At last, some may argue that the error may come from the coefficient n.This coefficient is empirical and we took n=0.25 because it gives the best accuracy in the results.And this value is identical for both winter wheat and maize datasets, which implies that the coefficient may have a universal value for all of the surfaces, but this still needs to be proved by more investigations.

Discussions
TSEBPS is proposed to estimate surface heat fluxes using TIR data obtained by space-borne sensors such as AVHRR, MODIS, etc.This kind of data is easily available and economical for the users, which is important for applications at regional or global scale with routinely schedule.For regional or global estimation of land surface evapotranspiration, sparsely vegetated surface is one of the situations of relatively larger uncertainty, where single layer model associated with TIR data can not simulate the canopy heat fluxes accurately.As a parameterization of the classical twolayer model, TSEBPS is reliable on the theory basis.It was shown in the evaluation using datasets over different vegetation canopies that TSEBPS-estimated evapotranspiration compared very well with the field measurement.The parameterization is based on the limiting cases of soil moisture, which is commonly accepted.The difference of TSEBPS is to consider foliage and soil independently at the limiting cases, and bring a key state of soil moisture into the model, i.e., transition state, which is based on the process of drying off after a rain or irrigation event when the soil surface is dry and the root zone is still wet.By the concept of transition state, we can hence define two different states of soil moisture in the canopy, i.e., before and after the transition, which represent the limit of soil moisture is only on Evaporation (E) or on both Evaporation (E) and Transpiration (T).The canopy heat fluxes are then easily predictable using the assumptions of the limiting cases associated with an interpolation method using TIR data.Commonly, all of the states of soil moisture can be described by such assumptions.However, there are exceptions when the soil surface is wet and root zone is relatively drier, which could be possible when there is heavy dew or light precipitation while the field has been under drought already.Under this circumstance, the relationship between surface temperature and soil moisture would be different from the assumption of TSEBPS, and TSEBPS-estimated heat fluxes would be of substantial error.Fortunately, this kind of exception is not a frequent event, i.e., once or twice during the whole growing season of crop, which will not affect the applicability of TSEBPS in the long term.
It can be found in the results that the TSEBPS-estimated heat fluxes under dense and wet canopy are similar to that under sparse and dry canopy with high accuracy.The empirical method that tries to relate TIR measurements with actual heat fluxes is able to produce good results and can be used widely for surfaces with different soil moisture and vegetation conditions.
According to the sensitivity analysis, the TSEBPSestimated heat fluxes are not sensitive to the assumed heat fluxes at the limiting cases as well as the error of the simulated temperature.On the other hand, it was found that higher accuracy can be obtained by using more complex model to allocate net radiation into soil and foliage.However, this could restrict the applicability of TSEBPS in satellite data.Compromise between accuracy and convenience has to be made.Fortunately, a simple method such as shown in Eqs. ( 4) to ( 6) can calculate soil and foliage net radiation reasonably and result acceptable heat fluxes in this study.On this meaning, the method proposed and used in this study is applicable for regional estimation of ET using satellite data.Results of evaluation of TSEBPS using satellite data will be reported by the authors in the near future.

Conclusions
Two-layer energy balance model has been validated and approved at many references.However, its application in remote sensing is still of problem because of short of component temperatures data.In this study, a parameterization scheme (TSEBPS) was proposed to utilize the two-layer model with traditional TIR observation data.The parameterization is based on the assumption of the changing process of sensible and latent heat fluxes of the foliage and sub-layer soil with the change of soil moisture at surface layer and root zone.The actual canopy heat fluxes are derived from the observed radiative surface temperature by comparing with the simulated temperatures at the limiting cases.Two datasets obtained in two different field experiments were used to evaluate the reliability of TSEBPS.The estimated canopy heat fluxes agreed well with the field measurements of heat fluxes.The uncertainties of the estimation are comparable to in-situ measurement uncertainties.The errors of TSEBPS mainly come from the following aspects, i.e., the assumption of the limiting cases, and the interpolation method of heat fluxes using the TIR observations.Although extensive evaluation should be carried out using more in-situ or remotely sensed data, the results of this study showed that the method proposed in this paper is reliable and can be used to estimate heat fluxes over sparsely vegetated surfaces.

Fig. 2 .
Fig. 2. Flow chart of the parameterization scheme of the two-layer models.

Fig. 4 .
Fig. 4. Time series of TSEBPS estimated latent heat flux bias (TSEBPS estimated minus measured latent heat flux) versus surface temperature gradient (radiative surface temperature minus air temperature) and wind speed.

Fig. 5 .
Fig. 5. Comparison between observations and TSEBPS modeled sensible and latentheat fluxes over maize canopy.Dashed line represents perfect agreement.

Table 1 .
Information about the turbulent and TIR measurements.

Table 2 .
Datasets used for the evaluation of TSEBPS.

Table 3 .
Statistics of TSEBPS estimated versus observed heat fluxes at site NW3 (RMSE: Root Mean Squared Error; MAD: Mean Absolute Deviation; MAPD: Mean Absolute Percentage Deviation; R 2 : coefficient of determination).

Table 4 .
Statistics of TSEBPS estimated versus observed heat fluxes at site NW4 (RMSE: Root Mean Squared Error; MAD: Mean Absolute Deviation; MAPD: Mean Absolute Percentage Deviation; R 2 : coefficient of determination).

Table 5 .
Statistics of TSEBPS estimated versus observed heat fluxes at site NW5 (RMSE: Root Mean Squared Error; MAD: Mean Absolute Deviation; MAPD: Mean Absolute Percentage Deviation; R 2 : coefficient of determination).

Table 6 .
Statistics of TSEBPS estimated versus observed heat fluxes at site YK (RMSE: Root Mean Squared Error; MAD: Mean Absolute Deviation; MAPD: Mean Absolute Percentage Deviation; R 2

Table 7 .
Statistics of TSEBPS estimated versus observed heat fluxes at three different growing stages of maize at site YK (RMSE: Root Mean Squared Error; MAD: Mean Absolute Deviation; MAPD: Mean Absolute Percentage Deviation; R 2 : coefficient of determination).