Parameterization of the Surface Energy Balance of a Shallow Water Table Grassland

: Extending instantaneous latent heat ﬂux to daily, monthly, or even yearly evapotranspiration (ET) is a fundamental issue in using remote sensing to estimate ET at local and regional scales. In this study, the extending parameterizations of the surface energy balance of a mid-latitude grassland with shallow water table (SWT) at diurnal and seasonal time scales are examined based on data measured by the eddy covariance system and automated weather station from Wageningen University from June 2014 to October 2018. The results show that the ratio of turbulent heat ﬂux to available surface energy (often called budget closure rate) ranges between 0.86 and 0.93 for warm times (March to October), and between 0.59 and 0.77 for cold times (November to February the following year). The parameterization models used to approximate the surface albedo and evaporative fraction (EF) are also evaluated. Although obvious variation under clear skies during daytime are observed, the constant EF and albedo method provided an acceptable estimation of the daily scale ET with an underestimation of about 6–8% for the grassland with SWT and parameterization of diurnal correction shows little improvement in both the bias and RMSE. The progression of daily ET shows a seasonal cycle, which follows the variation of the net radiation ﬂux. These results will be helpful for estimating ET at daily and long temporal scales based on satellite remote sensing.


Introduction
There is increased interest in evapotranspiration (ET) or consumptive use. ET is an important parameter in the process of surface water circulation and energy balance. Accurate estimation of surface ET is essential to better understand the global climate change, hydrometeorological processes, ecological environmental issues, agricultural irrigation, and watershed management [1]. ET can be estimated through field measurement using eddy covariance (EC), lysimeters, and the Bowen ratio [2][3][4]. The traditional field ET observations are usually limited in a few point locations and not sufficient to meet requirements of the regional scale applications. Fortunately, remote sensing provides a unique opportunity to retrieve a variety of surface parameters, such as reflectance and temperature that can be used to estimate the global and regional ET [5]. A lot of remote sensing models have been developed to predict ET, and they can be grouped into four main categories: empirical direct, residual, inference, and deterministic methods [6].
However, visible and infrared remote sensing only obtains the instantaneous observations when the satellite overpasses the surface at a certain time under clear conditions. In reality, the Earth's surface flux varies at different time scales. For the communities of hydrology, meteorology, and agriculture, the spatial and temporal continuous ET is more useful in practice [7]. Therefore, a fundamental problem in using remote sensing to estimate ET in local and regional scales is the extending of instantaneous latent heat flux to daily, monthly, or even yearly ET.
At present, there are many methods to do the temporal upscaling for the snapshot measurements, such as empirical model [8], sinusoidal relationship [9], evaporative fraction (EF) [10], reference EF [11], surface resistance [12], extraterrestrial radiation ratio [13], and so on. The most common one is the EF conservative method, which assumes that the daytime EF remains constant under the clear-sky condition, and then the total daily ET is estimated based on the available surface daily energy. The Surface Energy Balance Algorithm for Land (SEBAL) [14] and the Surface Energy Balance System (SEBS) [15] model achieve the upscaling from the instantaneous to daily scale by means of EF conservative base on the surface energy balance. However, there is still no consensus on the validity of the self-preservation of the EF.
Over the past few decades, many studies have examined the rationality of ET estimation under the constant assumption. Based on four days micrometeorological measurements during the First ISLSCP Field Experiment (FIFE), Shuttleworth et al. [10] found that the EF value changes slowly during daytime and the midday value can be used to estimate the daytime average. Crago [16] also used a few days of FIFE measurements to analyze the diurnal variation of EF. It shows that the EF has little variation during the daytime on clear days while fluctuating erratically on cloudy days. Van Niel et al. [17] reported that the EF self-preservation method leads to significantly underestimated of daily ET up to 34%. Peng et al. [18] used the EC system measurements from FLUXNET sites which cover different ecosystems and climates to examine the daytime EF constant hypothesis. They found that the EF between 11:00 and 14:00 LT (local time) agrees well with daytime (08:00 to 17:00 LT) EF on clear days while more variable on cloudy days. Leuning et al. [19] found that the EF varies with soil moisture and leaf area index during the day.
Apart from the EF, the surface albedo α, which affects the available energy, is another key variable to determine the ET estimation. The albedo also shows diurnal variation and it would lead to about 15% error when estimating the daily mean albedo from the 10:30 LT under the assumption that the albedo keeps a constant throughout the day [20].
In short, most of the above studies are usually based on measurements during a relatively short time period (no more than four seasons), which may not necessarily yield reliable results under different conditions. It remains to be a challenge to accurately depict how the varying EF and albedo influence the ET estimation. In an effort to address this question, the main objectives of this study included: (1) examining the patterns and parameterizations of the surface albedo and EF variation at hourly, daily, and seasonal scales for mid-latitude grassland with shallow water table (SWT); and (2) discussing the limitations and implications of this extending parameterizations approach. The ultimate objective in mind, which is not part of this work, is to use these extending parameterizations of the surface energy balance for estimating daily or long temporal based on satellite observations.

Physical Setting of Veenkampen Flux Site
Long-term observations were made at the Veenkampen Station by the Meteorology and Air Quality Group of Wageningen University from June 2014 to October 2018. The site (51.9814 • N and 5.6217 • E, about 5.0 m above sea level) is located in the center of Netherlands (see Figure 1). It consists of a grassland that is mainly covered with perennial species of ryegrass and rough bluegrass. During the growing season, from May to October, the grass cover is mowed weekly, thus it has a mean height of 0.1 m. The soil at the site is predominantly heavy basin clay resulting from the back-swamps of the Rhine River. There is a significant rainfall all year round with an average of 780 mm. Even the driest month still has more than 30 mm rainfall during the study period. It is noted that the site is about 300 m away from valley canal, which fits well with the hypothesis of an SWT. The SWT usually fluctuated between 0.2 m and 0.8 m under ground level during this study. However, it reaches a maximum depth of 1.2 m, which is far below the grass rooting depth during June 2018.
3 of 17 Figure 1. The Veenkampen Station equipped with the eddy covariance (EC) and automated weather station system, is located in Wageningen University at the center of Netherlands.

Field Measurements
The vertical turbulent heat fluxes (sensible (H) and latent heat (LE)) and agrometeorological parameters are measured by the EC system and automated weather station, respectively. The distribution of the instrumentation is shown in Figure 1 and the measurements are also listed in Table  1. The EC system was installed in a lattice tower at a height of 3.0 m. This system includes a threedimensional sonic anemometer (3-D Solent Res. Gill Instruments Ltd., model A1012R2, Lymington, United Kingdom), a fine wire thermocouple (home-made) and an open path infrared CO2 and H2O gas analyzer (IRGA) (LI-COR Inc., model LI-7500, Lincoln, NE, United States). The 3-D sonic anemometer and the IRGA are set 0.05 m apart. The raw data of the EC system sampled at 20.8 Hz were stored on a PC and processed at 30-min intervals with the EddyPro Software in advanced mode [21]. Quality flags were calculated for all the fluxes [22]. More details about the data processing can be found in Jacobs et al. [23]. The value "0" represents best quality fluxes, "1" represents questionable fluxes which suitable for general analysis such as annual budgets and "2" represents bad fluxes, which have been discarded from the raw dataset in this study.
Soil temperatures, Tso, are measured by Pt-100 elements at depths of 0.05, 0.10, 0.20, 0.50, 1.0 and 1.5 m; at 0.0625, 0.125, 0.25 and 0.50 m depths the soil moisture content was measured with a Time Domain Reflectometry (TDR) system (Campbell Sci., TDR 100, and probe type CS 610, Loughborough, UK). The soil heat flux is measured by a heat plate (TNO, WS 31-Cp, Delft, The Netherlands) buried at a depth of 0.06 m. The sunshine duration is measured by CSD3 Sunshine Duration Sensor (Kipp & Zonen, Delft, The Netherlands). The incoming (R ↓ ) and outgoing (R ↑ ) shortwave radiation fluxes are measured with an aspirated pyranometer (Kipp & Zonen, model CM11, Delft, The Netherlands) at a height of 1.5 m. The ratio of outgoing shortwave radiation to the incoming shortwave is used to calculate surface albedo. At the same height, the incoming (R ↓ ) and

Field Measurements
The vertical turbulent heat fluxes (sensible (H) and latent heat (LE)) and agrometeorological parameters are measured by the EC system and automated weather station, respectively. The distribution of the instrumentation is shown in Figure 1 and the measurements are also listed in Table 1. The EC system was installed in a lattice tower at a height of 3.0 m. This system includes a three-dimensional sonic anemometer (3-D Solent Res. Gill Instruments Ltd., model A1012R2, Lymington, United Kingdom), a fine wire thermocouple (home-made) and an open path infrared CO 2 and H 2 O gas analyzer (IRGA) (LI-COR Inc., model LI-7500, Lincoln, NE, United States). The 3-D sonic anemometer and the IRGA are set 0.05 m apart. The raw data of the EC system sampled at 20.8 Hz were stored on a PC and processed at 30-min intervals with the EddyPro Software in advanced mode [21]. Quality flags were calculated for all the fluxes [22]. More details about the data processing can be found in Jacobs et al. [23]. The value "0" represents best quality fluxes, "1" represents questionable fluxes which suitable for general analysis such as annual budgets and "2" represents bad fluxes, which have been discarded from the raw dataset in this study. The outgoing long wave radiation was used to evaluate the surface temperature by Stefan-Boltzmann law. The net radiation, R n is determined by summing up the total radiation budget. The long-wave radiometers, however, have a small sensitivity to shortwave radiation (15 W m −2 at full sunshine) and this thus affects the net radiation. A correction is therefore applied according to the specifications of the manufacturer (Kipp CG2 specifications) where R n * is the measured net radiation, α is the albedo of the underlying surface and R ↓max s is the maximum possible global irradiation, taken to be 1000 W m −2 .
An aspirated psychrometer measures the air temperature, Ta, and wet-bulb temperature, Tw, at a reference height zr = 1.5 m. At 0.10 m height, the air temperature Ta (0.10 m) is measured with a shielded Pt-100 thermometer (home-made). Note that the agrometeorological instruments are sampled at 0.25 Hz, and at 30-min intervals data are averaged.

Methods
Firstly, additional quality control is implemented to remove unreliable data. Rain periods and snowy days are discarded because of malfunction of the radiometers and EC system during these events. If the precipitation was recorded for more than one 60-minutes interval in a day or the albedo was larger than 0.4 and the daytime surface temperature less than 0 • C, that day is removed from further consideration. This exculpation is made empirically to avoid the dates in which the surface energy balance is complicated by precipitation (rainfall or snow). For each EC flux (LE and H), if more than 6 gaps on the 48 measurements during a day are presented, that date was removed from the set. Then the remaining gaps are filled by linear interpolation. The mismatch of measurement footprints between EC, radiation components, and potential satellite data is also checked. The integrated footprint was estimated after Schuepp et al. [24]. More than 95% of the cumulative flux stems from undisturbed grassed area (x < 690 m), except for the wind direction sector between 40 • and 80 • . From this sector, only 6% contributes to the flux and is ignored in this study.
Secondly, the surface energy budget (SEB) closure check is also performed to verify the quality of the data. The SEB closure is the base of remote sensing ET estimation. Although the storage terms are Water 2020, 12, 523 5 of 17 small in the mid-latitude grassland, they are not negligible for closing the SEB. The storage terms X can be estimated following Jacobs et al. [25]: where Sp (W m −2 ) is the photosynthesis flux, Sc (W m −2 ) is the crop enthalpy change, Sa (W m −2 ) is the air enthalpy change, Sd (W m −2 ) is the canopy dew water enthalpy change and Sq (W m −2 ) is the atmospheric moisture change. The storage changes are calculated over the measurement averaging time interval of 30 min. More details of the process can be found in Jacobs et al. [25]. According to the first law of thermodynamics, the sum of the estimated latent (LE) and sensible (H) heat flux is equivalent to all other energy sinks and sources where G0 is the surface soil heat flux (W m −2 ), Q the sum of all additional energy sources and sinks. G0 is approximated by a simple sinusoidal correction scheme. Typically, Q is neglected as a small term, and an imbalance between the remaining independently measured terms on the left-and right-hand sides of Equation (3) may indicate inaccurate estimates of the fluxes. The monthly SEB closure was evaluated using two different methods [26]. The first method is to derive linear regression coefficients from the least squares relationship between the half-hourly estimates of the dependent flux variables (LE + H) against the independently derived available energy (Rn − G0 − X); Another method to evaluate closure was to cumulatively sum Rn − G0 − X and LE + H over specified time periods, such as over an entire month or year during nocturnal periods, and calculate the energy balance ratio (EBR) The EBR gives an overall evaluation of energy balance closure at longer time scales by averaging over random errors in the half-hour measurements.
Thirdly, the temporal constancy of EF and albedo is verified during the daytime. EF and albedo are considered to be the most sensitive parameter in estimating flux. To be consistent with a number of previous studies that have been examined the variability of EF during the daytime of clear sky days, the analysis of EF and albedo are performed only for the period during 09:30 to 16:30 LT.
The total sunshine duration ratio (SDR, SDR = 1630 0930 sunshine duration record 7 hours ) is used to identify clear-sky days. In this study, only the days of which SDR equal to 1 are selected. There are 602 days under all-sky conditions of the final data set, and only 28 days under clear-sky conditions all day long. The instantaneous EF and albedo are calculated as where the subscript i refers to the i-th the instantaneous observation, R ↑ s i is the outgoing shortwave radiation and R ↓ s i is the incoming shortwave radiation. To investigate how well the instantaneous albedo to predict the daily mean albedo, a model proposed by Dickinson [27] is introduced: where α is the surface albedo, µ 0 is the cosine of the SZA, α(0.5) the albedo for µ 0 = 0.5 and d controls the strength of the SZA dependence. Four parameterizations of albedo are then compared with daily Thus Equation (5) can be rewritten as Fourthly, instantaneous derived results are compared against daily results. A number of sun-synchronous satellites, such as Landsat, SPOT, ERS − 1/2 and ENVISAT are widely used for ET estimation. They overpass the equator between 10:00 and 11:00 LT and the evaporative fraction calculated for 10:30 LT (EF 1030 ) is usually used to estimate the daily total latent heat flux in the current analysis. In this study, we are concerned about the daily evaporation, including the daytime and nighttime latent heat flux, although the nighttime latent heat flux is usually considered to be a minor component. To be useful for hydrology remote sensing, estimates of actual evaporation at the overpass time are compared with daily results, which includes both daytime and nighttime flux under all-sky conditions, not just clear-sky conditions during the daytime. The daily mean EF (EF 24h ) and albedo α d are defined as And then the daily total latent heat flux LE 24 where the EF i used to be considered as a constant during the day under the assumption that Rn and LE follow the same solar-related sinusoidal shape during the daytime cycle. To describe the increase of EF i during the afternoon, the latest version of SEBAL is correcting instantaneous EF i into a 24-h counterpart EF SEBAL i by means of an advection factor Ω that depends on the vapor pressure deficit (the difference between saturated vapor pressure, e 24 sat and actual vapor pressure, e 24 act ) following Hong et al. [28]: where Another parameterization method was proposed to correct the concave shape of EF during daytime to make a better estimation of daily ET by Hoedjes et al. [29].

Energy Balance Closure Check
Monthly energy balance is calculated using 30-minutes turbulent fluxes (LE + H) and available energy (Rn − G0 − X) components. Data are included in the analysis only when all the measurements are available. The nighttime (black dots) and daytime (grey dots) fluxes from all years combined for each month are shown in Figure 2. It shows that the nighttime observations at this station have no significant errors compared to observations in daytime, although the nighttime EC measurements are usually considered to be smaller and have relatively larger errors due to a stable atmospheric environment. The slopes of the regression through the origin range between 0.86 and 0.93 from March to October in warm times. It indicates that 86-93% of the variation in available energy (Rn − G0 − X) is accounted for the sum of the turbulent energy fluxes. These results are consistent with those reported in previous studies. It was reported that energy balance is in the range of 70-90% for a variety of ecosystems studied using EC measurements [30]. The sum of 30-minutes LE and H is strongly correlated with the sum of available energy (Rn − G0 − X) in these months, with a coefficient of determination R 2 value larger than 0.87, while the slopes decrease to a range of 0.59 and 0.77 for other months in cold times. Ideal closure is represented by an EBR of 1. The monthly EBR is ranging between 0.87 and 0.97 from February to October, which also closed to 1, while a spurious EBR much less or much greater than 1 appears from November to January of the following year in cold times.

Diurnal Variability of Surface Albedo
The surface albedo is another sensitive parameter in estimating flux. Figure 3 shows an example It shows that the nighttime observations at this station have no significant errors compared to observations in daytime, although the nighttime EC measurements are usually considered to be smaller and have relatively larger errors due to a stable atmospheric environment. The slopes of the regression through the origin range between 0.86 and 0.93 from March to October in warm times. It indicates that 86-93% of the variation in available energy (Rn − G0 − X) is accounted for the sum of the turbulent energy fluxes. These results are consistent with those reported in previous studies. It was reported that energy balance is in the range of 70-90% for a variety of ecosystems studied using EC measurements [30]. The sum of 30-minutes LE and H is strongly correlated with the sum of available energy (Rn − G0 − X) in these months, with a coefficient of determination R 2 value larger than 0.87, while the slopes decrease to a range of 0.59 and 0.77 for other months in cold times. Ideal closure is represented by an EBR of 1. The monthly EBR is ranging between 0.87 and 0.97 from February to October, which also closed to 1, while a spurious EBR much less or much greater than 1 appears from November to January of the following year in cold times.

Diurnal Variability of Surface Albedo
The surface albedo is another sensitive parameter in estimating flux. Figure 3 shows an example of 14 diurnal cycles of the incoming shortwave radiance (R ↓ s , black dots), the outgoing shortwave radiance (R ↑ s , blue dots) and the surface albedo (red dots) at the Veenkampen site during 2016 under clear sky conditions. For each day, the incoming and outgoing flux is nearly symmetric about noon, as expected. The albedo varies every day, ranging from 0.18 to 0.30. It always decreases with the local time in the morning, reaching its minimum near noon (~13:00 LT), and then increases again in the afternoon, as it is expected for a grass-covered surface from observations and the theory of Kimes et al. [31] that state, at smaller solar zenith angles (SZA), photons penetrate deeper into the vegetation and have multiple interactions with plant components and hence are more likely to be absorbed before escaping.  The Dickinson model (see Equation (7)) fits quite well for the 14 selected clear days during 2016 (green dots in Figure 3). Scatter plots for the 28 selected clear days during the whole study period between the observed daily mean clear-sky albedo and estimated daily mean albedo under those four assumptions are also shown in Figure 4. It shows that the constant-albedo assumption overestimates the daily mean albedo a variable amount that is close to 0.001 on a large proportion of days. The overall statistical results of the comparison between the estimated and observed daily mean albedo are list in Table 2. The estimation is biased positively because the albedo at 1030 LT is larger than time around noon, which radiation take a large proportion of the day. The empirical parameter d = 0.4 for grassland was presented by Briegleb et al. [32]. However, fixing d = 0.4 produced worse results than assuming constant albedo. Using the Dickinson model with d fitted to the data of the day, but setting the offset α(0.5) solely to match the 10:30 LT observation, the error in the estimate of daily mean albedo is lower. Assumption 4, using the Dickinson model with both parameters fitted to the day's albedo data, is the only case in which the model coefficients are adjusted to best reproduce the observed daily mean albedo. As expected, this gives the best result with errors being negligible. The mean difference between four estimated outgoing shortwave radiation (incoming shortwave radiation multiply by estimated albedo) and the measured outgoing shortwave radiation is −24.78 ± 13.85, −213.19 ± 81.11, −33.68 ± 14.43 and 2.85 ± 1.89 W m −2 , respectively. The Dickinson model (see Equation (7)) fits quite well for the 14 selected clear days during 2016 (green dots in Figure 3). Scatter plots for the 28 selected clear days during the whole study period between the observed daily mean clear-sky albedo and estimated daily mean albedo under those four assumptions are also shown in Figure 4. It shows that the constant-albedo assumption overestimates the daily mean albedo a variable amount that is close to 0.001 on a large proportion of days. The overall statistical results of the comparison between the estimated and observed daily mean albedo are list in Table 2. The estimation is biased positively because the albedo at 1030 LT is larger than time around noon, which radiation take a large proportion of the day. The empirical parameter d = 0.4 for grassland was presented by Briegleb et al. [32]. However, fixing d = 0.4 produced worse results than assuming constant albedo. Using the Dickinson model with d fitted to the data of the day, but setting the offset α(0.5) solely to match the 10:30 LT observation, the error in the estimate of daily mean albedo is lower. Assumption 4, using the Dickinson model with both parameters fitted to the day's albedo data, is the only case in which the model coefficients are adjusted to best reproduce the observed daily mean albedo. As expected, this gives the best result with errors being negligible. The mean difference between four estimated outgoing shortwave radiation (incoming shortwave radiation multiply by estimated albedo) and the measured outgoing shortwave radiation is −24.78 ± 13.85, −213.19 ± 81.11, −33.68 ± 14.43 and 2.85 ± 1.89 W m −2 , respectively. albedo is lower. Assumption 4, using the Dickinson model with both parameters fitted to the day's albedo data, is the only case in which the model coefficients are adjusted to best reproduce the observed daily mean albedo. As expected, this gives the best result with errors being negligible. The mean difference between four estimated outgoing shortwave radiation (incoming shortwave radiation multiply by estimated albedo) and the measured outgoing shortwave radiation is −24.78 ± 13.85, −213.19 ± 81.11, −33.68 ± 14.43 and 2.85 ± 1.89 W m −2 , respectively.

Evaporative Fraction
The The amplitude of LE is much larger than that of H, thus it seems less scattered than that of the LE in the same coordinate. The fluctuation is mainly due to the stochastic character of the turbulence process. To get a more stable EF, the instantaneous EF i is calculated by Equation (9) (13) and (16), respectively. The performance of the parameterization proposed by SEBAL and Hoedjes are compared with the observed and constant EF in Figure 5. It can be seen that assuming the EF self-preservation is not a valid hypothesis, since EF depicts a concave-up shape with a straight decrease in early morning and a sharp increase in the late afternoon. Because an increase in EF mainly results from an increase in incoming solar radiation and a decrease in atmospheric humidity, Hoedjes parameterizes the diurnal shape of EF as a function of incoming solar radiation and atmospheric humidity (see Equation (15)) [7]. On the contrary, the SEBAL EF parameterization is the only function of vapor pressure deficit (see Equation (14)), and for this reason, the Hoedjes parameterization approximates the observed EF diurnal variation in a better way than SEBAL.
The pattern of the seasonal variation of EF 24h is presented in Figure 6. Each of the points represents the daily average EF from 00:00 to 24:00 (Equation (10)). The EF 24h remains fairly constant at about 0.8 from 90 to 270 days of the year during the warm times, while it fluctuates dramatically around 0 to 1.5 during the cold times. The seasonal variation of EF 24h is a reflection of the climate of the area, particularly the solar radiation behavior. The fluctuations also appear in similar solar radiation conditions, with the rainfall and soil moisture as other possible reasons. Take the driest year 2018 in these years for example: it shows a little smaller EF 24h than others, especially in February, EF 24h less than 0.5 occurs.
The seasonal progression of EF 24h shows that, although the values are more stable in warm times than in cold times, there is still a day to day variation of approximately 0.3 in the EF in warm times in some occasions. The implication of this finding for the monitoring of EF is that it would be not sufficient to capture the evolution of EF using polar orbit satellites during warm times. Interpolation between the measurements cannot be done to estimate EF on days when no EF measurements are available. It means that, for remote sensing programs, the processing of daily images is necessary to estimate the variations of EF for large shallow water table mid-latitude grasslands in warm times. And proper temporal upscaling of the instantaneous satellite observations must be required to get a more accurate daily ET estimating [33]. compared with the observed and constant EF in Figure 5. It can be seen that assuming the EF selfpreservation is not a valid hypothesis, since EF depicts a concave-up shape with a straight decrease in early morning and a sharp increase in the late afternoon. Because an increase in EF mainly results from an increase in incoming solar radiation and a decrease in atmospheric humidity, Hoedjes parameterizes the diurnal shape of EF as a function of incoming solar radiation and atmospheric humidity (see Equation (15)) [7]. On the contrary, the SEBAL EF parameterization is the only function of vapor pressure deficit (see Equation (14)), and for this reason, the Hoedjes parameterization approximates the observed EF diurnal variation in a better way than SEBAL. The pattern of the seasonal variation of EF is presented in Figure 6. Each of the points represents the daily average EF from 00:00 to 24:00 (Equation (10)). The EF remains fairly constant at about 0.8 from 90 to 270 days of the year during the warm times, while it fluctuates dramatically  around 0 to 1.5 during the cold times. The seasonal variation of EF is a reflection of the climate of the area, particularly the solar radiation behavior. The fluctuations also appear in similar solar radiation conditions, with the rainfall and soil moisture as other possible reasons. Take the driest year 2018 in these years for example: it shows a little smaller EF than others, especially in February, EF less than 0.5 occurs. The seasonal progression of EF shows that, although the values are more stable in warm times than in cold times, there is still a day to day variation of approximately 0.3 in the EF in warm times in some occasions. The implication of this finding for the monitoring of EF is that it would be not sufficient to capture the evolution of EF using polar orbit satellites during warm times. Interpolation between the measurements cannot be done to estimate EF on days when no EF measurements are available. It means that, for remote sensing programs, the processing of daily images is necessary to estimate the variations of EF for large shallow water table mid-latitude grasslands in warm times. And proper temporal upscaling of the instantaneous satellite observations must be required to get a more accurate daily ET estimating [33]. The relationships between instantaneous (especially for the Landsat crossing time) and daily (both daytime and nighttime) EF under all weather conditions have been explored. All days (about 602 days(a)) separated into (b) warm times (90 to 270 DOY) and (c)cold times(the other days) are used, irrespective of weather conditions. Daily averaged EF (EF ) under all sky conditions are used to prepare scatter plots of the EF versus EF , depicted in Figure 7, to investigate the potential of using satellite remote sensing based data acquired during the morning hours. The 1:1 line (X = Y line) is shown as dashed. The blue solid line in each plot is the linear least squares regression fit to warm and cold times, respectively. There is not a strong relationship between EF and daily EF as expected. The coefficient of determination R 2 is about 0.6 for all days data. The R 2 for the regression lines are 0.53 and 0.59, while the root-mean-square error (RMSE) are 0.06 and 0.22 for warm and cold times, respectively. The 1:1 line shows that EF values are usually smaller than corresponding EF values during the warm times. EF larger than 1 occur in cold times, including January, February, November and December (see Figure 6). During these periods, when there is solar radiation deficit, the available energy is small, especially for the cloudy days. EF is expected to be lower at 1030 as compared with the rest of the day (Figure 7c). Meanwhile, EF at 1030 is usually larger than the EF when EF less than 1 occurs during cold times. This is because the magnitude of LE and H are very small and negative values maybe appear, which may lead to spurious EF. The relationships between instantaneous (especially for the Landsat crossing time) and daily (both daytime and nighttime) EF under all weather conditions have been explored. All days (about 602 days(a)) separated into (b) warm times (90 to 270 DOY) and (c)cold times(the other days) are used, irrespective of weather conditions. Daily averaged EF (EF 24h ) under all sky conditions are used to prepare scatter plots of the EF 1030 versus EF 24h , depicted in Figure 7, to investigate the potential of using satellite remote sensing based data acquired during the morning hours. The 1:1 line (X = Y line) is shown as dashed. The blue solid line in each plot is the linear least squares regression fit to warm and cold times, respectively. There is not a strong relationship between EF 1030 and daily EF 24h as expected. The coefficient of determination R 2 is about 0.6 for all days data. The R 2 for the regression lines are 0.53 and 0.59, while the root-mean-square error (RMSE) are 0.06 and 0.22 for warm and cold times, respectively. The 1:1 line shows that EF 1030 values are usually smaller than corresponding EF 24h values during the warm times.  Figure 8 displays the measured daily integrated ET flux plotted against the estimated results by (a) constant EF assumption set to 1030 LT and (b) corrected EF by Hoedjes method, Equation (16). It is striking that there are considerable scatters with an RMSE between observations and simulations of 0.65 MJ m −2 day −1 , and a determination coefficient R 2 of 0.95. Although it was expected that EF selfpreservation assumption would appear to lead to underestimates of ET due to the underestimate of EF as mentioned above, only for the larger ET shows obvious underestimate. Larger ET usually occurs in warm times (see Figure 9) when available energy is larger and constant EF of 10:30 LT tend to be smaller than most of the other EF during the day, leading to the obvious underestimated. In sum, the parameterization of EF corrected model accounts for the diurnal cycle and derive a more accurate estimate of daily ET (Figure 8b). The RMSE is reduced to 0.63 MJ m −2 day −1 and R 2 increased to 0.96. The slope of the linear fitted line is 1, which is an obvious improvement and perfect as compared to 0.92 when using constant EF. The daily mean latent heat flux measured from EC system and residual energy (RE, Equation (8)) is 5.10 and 5.59 MJ m −2 day −1 , respectively, as shown in Table 3. The constant EF and albedo method, which corresponds to the daytime overpass times of most polar-orbiting satellites at 1030 EF 24h larger than 1 occur in cold times, including January, February, November and December (see Figure 6). During these periods, when there is solar radiation deficit, the available energy is small, especially for the cloudy days. EF is expected to be lower at 1030 as compared with the rest of the day (Figure 7c). Meanwhile, EF at 1030 is usually larger than the EF 24h when EF 24h less than 1 occurs during cold times. This is because the magnitude of LE and H are very small and negative values maybe appear, which may lead to spurious EF.  (16). It is striking that there are considerable scatters with an RMSE between observations and simulations of 0.65 MJ m −2 day −1 , and a determination coefficient R 2 of 0.95. Although it was expected that EF self-preservation assumption would appear to lead to underestimates of ET due to the underestimate of EF as mentioned above, only for the larger ET shows obvious underestimate. Larger ET usually occurs in warm times (see Figure 9) when available energy is larger and constant EF of 10:30 LT tend to be smaller than most of the other EF during the day, leading to the obvious underestimated. In sum, the parameterization of EF corrected model accounts for the diurnal cycle and derive a more accurate estimate of daily ET (Figure 8b). The RMSE is reduced to 0.63 MJ m −2 day −1 and R 2 increased to 0.96. The slope of the linear fitted line is 1, which is an obvious improvement and perfect as compared to 0.92 when using constant EF.

Evapotranspiration Flux
The daily mean latent heat flux measured from EC system and residual energy (RE, Equation (8)) is 5.10 and 5.59 MJ m −2 day −1 , respectively, as shown in Table 3. The constant EF and albedo method, which corresponds to the daytime overpass times of most polar-orbiting satellites at 1030 LT, is used to estimate the daily ET. Little underestimation about 6-8% of the daily ET is generally found in both the EC and RE methods. The absolute bias decreased by 3-6 % when the parameterization of diurnal correction was applied. The RMSE reduced by 0.09 MJ m −2 day −1 and 0.02 MJ m −2 day −1 for EC and RE methods with diurnal corrected. It could be concluded that the conventional constant EF and albedo method provided an acceptable estimation of the daily scale for the shallow water table grassland and the parameterization of diurnal correction shows little improvement in both the bias and RMSE. occurs in warm times (see Figure 9) when available energy is larger and constant EF of 10:30 LT tend to be smaller than most of the other EF during the day, leading to the obvious underestimated. In sum, the parameterization of EF corrected model accounts for the diurnal cycle and derive a more accurate estimate of daily ET (Figure 8b). The RMSE is reduced to 0.63 MJ m −2 day −1 and R 2 increased to 0.96. The slope of the linear fitted line is 1, which is an obvious improvement and perfect as compared to 0.92 when using constant EF. The daily mean latent heat flux measured from EC system and residual energy (RE, Equation (8)) is 5.10 and 5.59 MJ m −2 day −1 , respectively, as shown in Table 3. The constant EF and albedo method, which corresponds to the daytime overpass times of most polar-orbiting satellites at 1030   The progression of daily ET, applying the Dickinson albedo and Hoedjes EF parameterizations, and corresponding measured net radiation flux and water level during June 2014 to October 2018 are presented in Figure 9. It shows a general trend of daily ET dominated by a seasonal cycle, which follows the variation of net radiation flux. The daily ET increases from January to April, reaching its peak during May, June, and July, and then decreasing from August to December. The daily ET is usually less than 1 mm/day during the cold times (DOY 1 to DOY 90 and DOY 271 to DOY 365), thus lead to the cumulative ET in the warm times (DOY 91 to DOY 270) that takes about 88% of the annual value. It is noted that there are many clouds in June and July, which lead to the fluctuations of net radiation flux consequent to the ET, although the extraterrestrial radiance is the largest of the year. The water table is always shallow and fluctuated between 0.2 m and 0.5 m in cold times and between 0.3 m and 0.8 m in warm times. However, for a severe drought month during this study, June 2018, a maximum depth of 1.2 m is reached, which is far below the grass rooting depth. Thus even if the net radiation flux was sufficient, the ET was smaller than that in the same days of other years.

Limitation of This Approach
The SEB closure is the base of remote sensing ET estimation. The G0 is considered to be small in magnitude and a simple sinusoidal correction scheme is applied to approximate the soil heat flux at the surface in this study. Jacobs et al. [25] improved the surface energy budget closure rate by 9% using the harmonic analysis technique. It would be a promising way to improve our SEB closure. This study mainly investigated the parameterizations of EF and albedo to account for the diurnal and seasonal variation of ET estimation. However, it does not account for the direct stomatal response to water availability. Stomatal regulation depends on volumetric water content (VWC). According to the Soil & Water Assessment Tool (SVAT) simulations [34], the ET is water-limited, leading to a decorrelation with solar incoming radiation when VWC is less than 0.2 (m 3 m −3 ) for the grass cover, which had a leaf area index of about 2.9 at this site. The physiological control in the grassland is mild in potential conditions, or the net effect of EF increase due to lower RH values and stomatal closure in the afternoon. Although the shallow water table and continuous intermittent rainfall at this site usually leads to a high VWC (larger than 0.3 m 3 m −3 ), VWC is consistently close to 0.2 m 3 m −3 for several weeks during the extremely dry seasons in the Netherlands in 2018 (see Figure 10). Therefore, parameterization of stomatal regulation also needed to be considered for this site under the low VWC conditions in future works.
It is considered that the surface albedo only depends on the SZA and free of the azimuth under clear-sky conditions. So, the albedo should be equal at the same SZA before and after noon. However, the observations in Figure 3 indicate that there is an asymmetry in the morning and afternoon, and generally the morning albedo is little smaller than the afternoon one. The diurnal variability of dew, frost, evaporation, and wind may be the cause of the asymmetry [35]. For the selected data after 10:30 LT during warm times, the dew and frost may not exist, thus the wind may be the possible reason. At Veenkampen, for the selected days, the prevailing wind was eastern, which could cause the grass to lean toward the west and shift the albedo minimum into the morning. generally the morning albedo is little smaller than the afternoon one. The diurnal variability of dew, frost, evaporation, and wind may be the cause of the asymmetry [35]. For the selected data after 10:30 LT during warm times, the dew and frost may not exist, thus the wind may be the possible reason. At Veenkampen, for the selected days, the prevailing wind was eastern, which could cause the grass to lean toward the west and shift the albedo minimum into the morning.

How to Handle under Cloudy Conditions
The high temporal resolution of geostationary satellites provides a unique opportunity to handle the ET estimation under cloudy conditions. Although the instantaneous geostationary satellites observation is also affected by the cloud, with the movement of the cloud, it can effectively improve the availability of the observations and overcome the influence of the cloud. One or more instantaneous observations of geostationary satellites are then combined into a parametric model to give an accurate estimate of the day or longer time scale ET. However, the kilometric spatial resolutions severely limits water management at the field scale. Multi-sensor data fusion would be a promising way to solve the problem. Spatial-Temporal Adaptive Reflectance Fusion Model (STARFM) [36] provides a model to fuse MODIS (coarse spatial resolution, 250 m to 1 km) and Landsat (medium spatial resolution, 30 m) data. With the rapid development of remote sensor and computer technology, the spatial resolution and accuracy of geostationary satellites are also increasing, which is promising to apply the fusion model, such as STARFM, to get a high spatial and temporal resolution data.

Conclusions
The assumption of constant EF and albedo are commonly used to extrapolate the daily values from the instantaneous observation by remote sensing. However, evidence for this constant assumption method is based on limited duration field measurements. Taking advantage of long-term ground-based measurements from Veenkampen, which locates at the mid-latitude grassland with SWT in Netherland, the parameterization of the surface energy balance is examined in this study. The following are the main conclusions that can be drawn from the study: