**4. Results**

#### *4.1. Land Surface Heterogeneity at FluxNet Sites*

Figure 4a shows an 18 km subset TM albedo for IT-Col and a magnified image of the measurement point and spatial boundaries of 1.0 km<sup>2</sup> and 1.5 km2, which provides a detailed visual representation of the landscape heterogeneity at the site. Figure 4b shows the variogram estimator (point values) for the TM albedo subsets. The spatial variability is more obvious over the larger squared regions (1.5 km<sup>2</sup> and 2.0 km2), and different land cover types are more likely to be found at larger separation distances. This variogram estimator reaches an asymptote or constant variance between spatial uncorrelated samples near the sample variance. We calculated the range and Rcv for each Flux site to determine the heterogeneity of MODIS pixels as summarized in Table 2.

When the TM footprint increased from 1.0 km to 1.5 km, Rcv increased 0%–10% over CA-Oas, IT-Col, IT-Ro2, FR-Pue, AU-Wac, CH-Dav, NL-Loo, DE-Kli, IT-BCi, FR-Gri, US-ARM, AU-DaP, and US-IB2. These sites differ less from the surrounding landscape; when the TM footprint changes, Rcv only changes slightly relative to the smaller landscapes in the surrounding regions. The size of Rcv depends on the degree of landscape similarity within the given range. When surface heterogeneity around the site is weak, the overall degree of stationarity between regions is similar, the overall variability is small, and the value of Rcv is relatively low. Rcv increased 32.46% over HI-Fyy and 86.38% over IT-Sro as the TM footprint increased from 1.0 km to 1.5 km. Based on the range of FI-Hyy beyond the 1 km<sup>2</sup> limit, the landscapes around the two sites are obviously different; the small lakes near FI-Hyy are not within the 1 km<sup>2</sup> boundary, IT-Sro is near the sea, and internal (CV1 km) regions do not include the sea. At the US-Bar, MY-Pso, and AU-Stp stations, the results for RCV were moderate (10%–20%). There are two types of vegetation around the three sites which include mixtures of grass and trees. There are some grasslands located approximately 1.0 km north of the US-Bar tower and southeast of the MY-Pso tower. A broadleaf forest is located southwest of the AU-Stp station, but the internal (CV1 km) does not include the forest. The overall fluctuation between regions was not high—the three sites were moderately heterogeneous.

**Figure 4.** *Cont*.

**Figure 4.** (**a**) Subset of 18.0 km TM albedo at IT-Col (measurement point) and spatial boundaries of 1.0 km<sup>2</sup> and 1.5 km2. (**b**) Variogram plot, spherical model (dotted curves), and sample variance obtained via TM albedo on 25 July, 2009 for 1.0 km<sup>2</sup> (asterisks), 1.5 km<sup>2</sup> (diamonds), and 2.0 km<sup>2</sup> (squares) regions.


**Table 2.** Rcv for selected ground sites.

#### *4.2. Estimating and Verifying Single-Point Time Series Albedo*

We included five vegetated land surface types in this study: Deciduous broadleaf forest, evergreen broad-leaved forest, evergreen coniferous forest, grassland, and farmland. Three sites were subjected to validation for each land surface type. The results are shown in Figure 5.

The algorithm starts with cloudless TM albedo data spanning the whole year. The starting days for each site are different, so the MODIS albedo data for different land surface types have different degrees of overestimation and underestimation; however, they are relatively continuous in time and all reflect the variation characteristics of surface albedo in snow-free periods. The EnKF algorithm combines the advantages of MODIS and TM sensors. The estimated results were very close to the directly calculated TM albedo data when adequate TM data were included. When TM data is absent, the albedo background is used as the final estimation to maintain a complete time series. Our results also showed that with the continuous introduction of TM observations, the errors induced by MODIS background data could be corrected by EnKF. The albedo background is generated by resampling MODIS albedo data.

We also evaluated the heterogeneity of the test sites (Section 4.1) within the entire MODIS 500 m spatial scale. The results indicated that the proposed algorithm can adapt to heterogeneous surfaces (Figure 5e,h,o).

Figure 6 shows scatter plots of the estimated albedos and ground measurements. Different colors represent the three different sites for each land surface type.

As shown in Figure 6, RMSEs of the estimated albedos in deciduous broadleaf forest (DBF), evergreen broadleaf forest (EBF), evergreen needleleaf forest (ENF), cropland (CRO), and grassland (GRA) were 0.0152, 0.0085, 0.0087, 0.0152, and 0.0109; the R<sup>2</sup> values were 0.5980, 0.7061, 0.9190, 0.8191, and 0.6911, respectively. The RMSEs for all land surface types were less than 0.02, which meets the requirements of global and regional climate models [53]. The correlation coefficients between the estimated surface albedos and the surface measurements were high ( R<sup>2</sup> > 0.6). The estimation results were close to the line *y* = *x*. The accuracies of ENF and CRO were the highest among all land surface types with R<sup>2</sup> values of 0.9190 and 0.8191, respectively. The lowest R<sup>2</sup> value was 0.5980 for DBF, which can be attributed to the dispersion of TM data (Figure 1). At ENF and CRO sites, the distribution of TM data was more uniform across the whole time series. At the DBF site, there was a slightly longer period of missing data which prevented the albedo from being updated in a timely manner. The estimation accuracy depended on the temporal distribution of the available TM data at a certain degree, as temporally uniform distributed TM data led to higher estimation accuracy. Nevertheless, the proposed algorithm greatly improves the estimation accuracy compared to field observation alone.

**Figure 5.** *Cont*.

**Figure 5.** Time series albedo estimations based on EnKF for DBF (**<sup>a</sup>**–**<sup>c</sup>**), EBF (**d**–**f**), ENF (**g**–**i**), CRO (**j**–**l**), GRA (**<sup>m</sup>**–**<sup>o</sup>**); three sites for each surface type. Grey vertical line represents error between estimated and measured values.

**Figure 6.** Estimated albedos and ground measurements at AmeriFlux and Flxunet sites: (**a**) DBF; (**b**) EBF; (**c**) ENF; (**d**) CRO; (**e**) GRA. Different colors represent different sites.

#### *4.3. Regional Timing Albedo Estimation and Verification*

We validated the proposed assimilation algorithm on five types of land surface areas, each area containing a ground-based station. The five stations, as mentioned above, were US-Bar, AU-Wac, IT-Sro, US-Arm, and AU-Stp (Table 1). The results are shown in Figure 7. The daily estimation result was difficult to display and the beginning date for each site was different, so we used nine days-worth of estimation results from each site.

The EnKF assimilation method uses TM albedo as observation data. The estimation is strongly dependent on the quantity of TM albedo values. To obtain time-continuous, high spatial resolution albedo data, the EnKF algorithm updates the background field from the first cloudless TM image; therefore, assimilation results will be improved as more TM images become available. Figure 7 shows the estimated time series albedos for five different land cover types. Again, the starting point of each time series was not exactly the same as was determined by the availability of high-quality TM data. Starting with the first cloudless TM image, the algorithm can obtain an albedo dataset with a temporal resolution of one day and a spatial resolution of 30 m (shown here as 25-day or 35-day intervals), thereby reflecting albedo fluctuations within a one-year period.

As verification, we compared the albedo generated by the EnKF algorithm and the TM albedo directly estimated on the same date. Figure 8 shows that the EnKF and TM albedo results are consistent. The scatter plot has RMSE values between 0.0031 and 0.0112, and R<sup>2</sup> values between 0.8584 and 0.9494. One exception was the farmland area (US-Arm site), where regional estimation and TM albedo verification results had R<sup>2</sup> lower than the other four surface types, and RMSE higher than other landforms. This may be due to the fact that farmland surface vegetation is more affected by anthropogenic factors than other land surface types; the jump in surface albedo makes the deviation of farmland albedo estimations larger than others. In this case, inducing more high-resolution images into the assimilation process could effectively improve the albedo estimation accuracy.

**Figure 8.** Shortwave albedo map generated by EnKF method with directly estimated TM albedo for five different land surface types: (**a**) DBF at US-Bar site (2009); (**b**) EBF at AU-Wac site (2007); (**c**) ENF at IT-Sro site (2009); (**d**) CRO at US-Arm site (2010); (**e**) GRA at AU-Stp site (2009). Right-hand scatter plots show consistency between the two sets of results, as number of dots gradually decreases from red to blue.

The DBF, EBF, and CRO sites selected in this study are heterogeneous research areas and have RMSE values less than 0.01 and R<sup>2</sup> values greater than 0.9, which indicate that the algorithm is well-suited to albedo estimation in heterogeneous areas.
