# Surface albedo and toc-r 300 m products from PROBA-V instrument in the framework of Copernicus Global Land Service

Jean-Louis Roujean, Jonathan Leon-Tavares, Bruno Smets, Patrick Claes, Fernando Camacho de Coca, Jorge Sanchez-Zapero

► To cite this version:

Jean-Louis Roujean, Jonathan Leon-Tavares, Bruno Smets, Patrick Claes, Fernando Camacho de Coca, et al. Surface albedo and toc-r 300 m products from PROBA-V instrument in the framework of Copernicus Global Land Service. Remote Sensing of Environment, 2018, 215, pp.57-73. 10.1016/j.rse.2018.05.015. hal-02354384

HAL Id: hal-02354384

https://hal.science/hal-02354384

Submitted on 28 Aug 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Surface albedo and toc-r 300 m products from PROBA-V instrument in the framework of Copernicus Global Land Service

<figure><img src="image_3.png" /></figure>

Check for
updates

Jean-Louis Roujean<sup>a,*</sup>, Jonathan Leon-Tavares<sup>b</sup>, Bruno Smets<sup>b</sup>, Patrick Claes<sup>b</sup>,
Fernando Camacho De Coca<sup>c</sup>, Jorge Sanchez-Zapero<sup>c</sup>

## ARTICLE INFO

### Keywords:

Albedo

Vegetation

Satellite

Time series

## ABSTRACT

PROBA-V instrument launched in 2013 is offering a global daily coverage at pixel resolutions of 333 m and 1 km in three spectral bands (BLUE, RED, NIR) and 600 m for shortwave infrared (SWIR). The PROBA-V mission is the follow-on of the VEGETATION program started in 2000, which allowed generating long-term series at 1 km pixel resolution. The PROBA-V products belong to the Copernicus Global Land Service portfolio (http://land.copernicus.eu/global/). The sensor design of PROBA-V with oriented cameras offers a wide field of view (FOV) for sampling the BRDF (Bidirectional Reflectance Distribution Function). This paper details the methodology implemented at the premises of VITO (Flemish Institute for Technological Research) with the aim to disseminate routinely from PROBA-V daily observations for both surface albedo (SA) and top-of-canopy corrected reflectance (TOC-R) products. The method classically operates a selection of cloudless scenes, performs atmospheric corrections, and finally applies a correction of directional effects on a pixel per pixel basis. The synthesis period is the decade and the composite period is 20 days. Such choice is a pointwise sampling as being a trade-off between the availability of clear scenes and the timescale for phenology. Regarding the albedo catalogue, a narrow-band to broadband conversion is stipulated. A recurrent technique serves for gap-filling based on the spread of weighed a priori data. Additional information concerns the quality flag and the age of the product. Preliminary accuracy assessment is performed through a comparison with the Moderate Imaging Spectroradiometer (MODIS) Collection 6. Dependable spatial consistency is reached except for wintertime with deviations in terms of rmse (root mean square errors) about 0.03 for visible and shortwave domains, and 0.04 for near infrared. Besides, both PROBA-V and MODIS C6 exhibit close time profiles, marked by smoothness or rapid transitions. Results over 10 confidence sites reveals rmse values of 0.032 and bias of 0.01 over the 2014 full annual cycle.

## 1. Introduction

Land surface albedo is the cornerstone for characterizing the energy balance in the coupled surface-atmosphere system and also constitutes an indispensable input quantity for soil-vegetation-atmosphere transfer models. It yields an Essential Climate Variable (ECV) as established by the Global Climate Observing System (GCOS) (GCOS, 2016) with given guidelines for its long-term validation (http://www.qa4ecv.eu/ecv/albedo). Knowing the surface albedo, the net radiation at the surface can be estimated and besides the whole energy budget. Heretofore, three spectral broadband ranges, namely the solar spectrum (400–3000 nm), the visible (400–700 nm) and the near- and shortwave-infrared (700–3000 nm), were deemed the relevant quantities. Actually,

any change in the short-wave (solar) albedo can be tenuous because of the counter-balancing between broadband visible and near infrared surface albedo. This fully justifies the dissemination of the three broadband albedo products although one could be derived from the two others. Noteworthy, the spectral range for visible broadband is matching with PAR (Photosynthetically Active Radiation) range to depict the carbon budget. As vegetation absorbs most of the PAR radiation, therefore PAR albedo is particularly sensitive to greenness. On the other hand, near-infrared albedo is high for leafy vegetation and low for woody material comparatively to visible albedo.

Actually, there exists variant definition of albedo products according to the domain of directional integration (Schaepman-Strub et al., 2006). It notably places a regard to the fraction of direct versus

<footnote><sup>a</sup> UMR3589-CNRM, CNRS/METEO-FRANCE, 42 avenue Coriolis, 31057 Toulouse, France</footnote>

<footnote><sup>b</sup> VITO, Boeretang 200, 2400 Mol, Belgium</footnote>

<footnote><sup>c</sup> EOLAB, Valencia, Spain</footnote>

<footnote>* Corresponding author at: UMR-5126, CESBIO, 18 avenue E. Belin, bpi 2801, 31401 Toulouse, France.
E-mail address: jean-louis.roujean@cesbio.cnes.fr (J.-L. Roujean).</footnote>
diffuse solar radiation. The salient albedo products are the Directional-Hemispherical Reflectance (DHR) – also called Black Sky Albedo (BSA) - and the Bi-Hemispherical Reflectance (BHR) - also called White Sky Albedo (WSA). Their combination relative to the ratio of sky irradiance leads to the so-called Blue Sky Albedo. Actually, the Blue-Sky Albedo is the true albedo to be measured in situ.

The resolution of 333 m offered by PROBA-V sensor will prompt new applications in domains encompassing agriculture, forestry, land use, land cover, hydrology and weather forecasts areas. The point-wise spatial-temporal resolutions of PROBA-V is also prone to leverage geo-engineering activities in order to dampen the effects of a changing climate. The Copernicus Global Land Service (CGLS) operates “a multi-purpose service component” that offers a series of bio-geophysical products on the status and evolution of land surface at global scale (http://land.copernicus.eu/global). Timely production and delivery of set of parameters exacerbate the constitution of long term series of satellite-based products elaborated in a coherent manner. The primary objective of CGLS is to continuously monitor the status of land territories and to supply reliable geo-information to decision makers, businesses and citizens to define environmental policies and take right actions. ImagineS (Implementing Multi-scale Agricultural Indicators Exploiting Sentinels) project from FP7 (Framework Program Seventh) was at the root of the development of cutting-edge retrieval methods of key biophysical variables, amongst which the land surface albedo. The algorithm to measure the land surface albedo from PROBA-V is a trimmed methodology previously implemented in operational for Meteoat Second Generation (MSG) (e.g. Geiger et al., 2008). The insurance of the continuity with past product from SPOT/VEGETATION at 1 km is enacted by the follow-on dissemination of 1 km PROBA-V product, owing to CGLS. The algorithm has been fine-tuned with time for the purpose of an enhanced efficient computation and service requirements.

The paper first reviews the background theory about the kernel-
driven BRDF (Bidirectional Reflectance Distribution Function) ap-
proach. The BRDF model parameters serve to estimate both spectral
albedos and TOC-R. A narrow- to broadband conversion is then per-
formed. If most satellite projects adopted the kernel-based approach, a
variant was operated for MISR (Multi-angle Imaging Spectro-
Radiometer), further assessed in terms of noise (Lucht and Lewis,
2000), then of added-value for MODIS (Jin et al., 2002) and also in
virtue of its potential to map snow albedo (Stroeve and Nolin, 2002).
Whilst a surface albedo product is useful for surface energy balance and
radiation forcing at surface level, TOC reflectance normalized to geo-
metry of reference is more inclined to serve for the monitoring of the
surface resources and the derivation of vegetation indices. Section 2
presents the instrument, the calibration accuracy, the levels of data
processing and quality. Section 3 tells about the methodology im-
plemented. Section 4 details the tools and criteria for products eva-
luation. Section 5 provides a preliminary assessment of the products
quality. Section 6 concludes the study and stresses future prospects.

## 2. Characteristics of PROBA-V instrument

### 2.1. *Principle of measurement*

PROBA-V payload named VGT was launched in 2013 for 7 years and is fully comparable to the previously VEGETATION sensor embarked on SPOT (Satellite Probatoire d'Observation de la Terre). It is a multi-spectral push-broom spectrometer. The payload consists of three identical cameras, equipped with a very compact Three Mirror Anastigmat (TMA) telescope. Each TMA has a FOV of 34° with four spectral bands. The limit of view zenith angle is 75° (e.g. http://proba-v.vgt.vito.be/sites/proba-v.vgt.vito.be/files/Product_User_Manual.pdf). Three spectral bands belong to the visible range (460 nm for Blue, 658 nm for Red and 834 nm for NIR) plus a SWIR band (1610 nm). VGT is restricted to imaging land and dedicated calibration zones. Each camera owes its own land sea mask that allows removing the sea pixels. About 14 near

polar orbits per day are registered. PROBA-V flies at 820 km altitude. The swath width of 2250 km ensures a daily coverage of land masses above 35° latitude with however a limitation to 75° North and 56° South. About 90% daily coverage is obtained in the equatorial zones. The Ground Sampling Distance is 100 m (VNIR) and 200 m (SWIR) at nadir and 360 m (VNIR) and 690 m (SWIR) at the edge of the swath.

## 2.2. *Data preprocessing and performances*

The images are projected in the grid plate carrée for Level 1-b and the geodetic datum is WGS84. The pixel co-ordinates are given for the center of the pixel. For details on the radiometric performances of PROBA-V, we will refer the reader to the dedicated link (https://earth.esa.int/web/sppa/mission-performance/esa-3rd-party-missions/proba-v/products-and-algorithms/products-information). To be outlined here that the PROBA-V S1 Top of Canopy (TOC) reflectance values are synthesis of the pixels from the three cameras having harmonized spectral responses functions. Daily PROBA-V composites (S1) of TOC reflectance, at a spatial resolution of 333 m, are the primary sensor data serving as input for both TOC-R and surface albedo algorithms. Note that official references to PROBA-V products indicate 300 m although the true resolution is 333 m. In case of multiple observations per day, the maximum of Normalized Difference Vegetation Index (NDVI) enacts a criterion of selection. The current PROBA-V cloud detection method (implemented in collection 1, C1) identifies the presence of clouds based on land cover class, climatology background surface reflectance per pixel and a set of decision rules (Sterckx et al. (2014), Dierckx et al., 2014, Wolters et al., 2017). The cloud mask is majorly inherited from VEGETATION. The Digital Elevation Model is from GTOPO30 (U.S. Geological Survey). The columnar water vapor and ozone contents, also the atmospheric pressure, are input fields issued from the numerical weather prediction model of the European Centre for Medium-Range Weather Forecasts (ECMWF). The ground segment of PROBA-V applies an atmospheric correction on a pixel-per-pixel basis for cloud-free pixels. The atmospheric correction is performed using the SMAC (Simplified Method for Atmospheric Correction) software (Rahman and Dedieu, 1994). SMAC considers linear parameterizations for the absorption and scattering components of molecules and aerosols based on the physics of the 6S code. These relationships include tunable coefficients depending on aerosol type and PROBA-V channels. However, a continental type is taken everywhere whereas the aerosol optical depth (AOD) at the wavelength of 550 nm is prescribed as a function of latitude. The S1 TOC PROBA-V reflectance is distributed with the Status Map (SM), which tells about the quality of the product (radiometry quality, cloud mask, etc.) (see Table 1). A land sea mask (LSM) is used to delineate the coastline. Note that for time being, the inner water bodies are not masked by the application of LSM.

## 3. Methodology description

### 3.1. *Algorithm overview*

The operational processing scheme of the land surface albedo and
TOC-R algorithm is depicted in the flow chart of Fig. 1. It encompasses
three successive steps: the spectral TOC reflectance values serve as the
input quantities for the inversion of a linear kernel-driven BRDF model,
which allows taking into account the angular dependence of the re-
flectance factor. A well-established approach for an operational com-
putation of the surface albedo is based on semi-empirical BRDF kernel
model. Such category of models has received a great deal of attention
and effort from the optical remote sensing community in the last dec-
ades (Roujean et al., 1992; Barnsley et al., 1994; Wanner et al., 1995;
Strahler, 1994; Hu et al., 1997). The approach is based on a decom-
position of the bi-directional reflectance factor into a number of kernel
functions which are associated to the dominant light scattering pro-
cesses, e.g. geometric and volumetric effects, a separation between the
<table><caption>**Table 1**
Explanation of the pixel quality indicators in the Status Map Dataset. Bits indicated with an asterisk are only available for Level2A data.</caption><thead><tr><th>Bit</th><th>Description</th><th>Value</th><th>Key</th></tr></thead><tbody><tr><td rowspan="5">Bits 0–2</td><td rowspan="5">Cloud/Ice Snow/Shadow Flag</td><td>000</td><td>Clear</td></tr><tr><td>001</td><td>Shadow</td></tr><tr><td>010</td><td>Undefined</td></tr><tr><td>011</td><td>Cloud</td></tr><tr><td>100</td><td>Ice</td></tr><tr><td rowspan="2">Bit 3</td><td rowspan="2">Land/Sea</td><td>0</td><td>Sea</td></tr><tr><td>1</td><td>Land</td></tr><tr><td rowspan="2">Bit 4</td><td rowspan="2">Radiometry quality SWIR flag</td><td>0</td><td>Bad</td></tr><tr><td>1</td><td>Good</td></tr><tr><td rowspan="2">Bit 5</td><td rowspan="2">Radiometry quality NIR flag</td><td>0</td><td>Bad</td></tr><tr><td>1</td><td>Good</td></tr><tr><td rowspan="2">Bit 6</td><td rowspan="2">Radiometry quality RED flag</td><td>0</td><td>Bad</td></tr><tr><td>1</td><td>Good</td></tr><tr><td rowspan="2">Bit 7</td><td rowspan="2">Radiometry quality BLUE flag</td><td>0</td><td>Bad</td></tr><tr><td>1</td><td>Good</td></tr><tr><td rowspan="2">Bit 8*</td><td rowspan="2">SWIR coverage</td><td>0</td><td>No</td></tr><tr><td>1</td><td>Yes</td></tr><tr><td rowspan="2">Bit 9*</td><td rowspan="2">NIR coverage</td><td>0</td><td>No</td></tr><tr><td>1</td><td>Yes</td></tr><tr><td rowspan="2">Bit 10*</td><td rowspan="2">RED coverage</td><td>0</td><td>No</td></tr><tr><td>1</td><td>Yes</td></tr><tr><td rowspan="2">Bit 11*</td><td rowspan="2">BLUE coverage</td><td>0</td><td>No</td></tr><tr><td>1</td><td>Yes</td></tr></tbody></table>

<figure><img src="image_4.png" /><figcaption>**Fig. 1.** Flow chart of the algorithm for BRDF model inversion and albedo and TOC-R determination.</figcaption></figure>

soil and vegetation, or the conjunction between media which are optically thick and thin (Lucht and Roujean, 2000).

The kernel-based approach was adopted to map the surface albedo for a number of space-borne sensors like POLDER, SeaWiFS, VEGETATION, MODIS and MERIS (e.g., Leroy et al., 1997; Justice et al., 1998; Wanner et al., 1997; Strahler et al., 1999; Samain et al., 2006; Baret et al., 2007; Muller et al., 2012). The Top-Of-Canopy Reflectance (TOC-R) product is however not always an output product. A TOC-R product is deemed relevant to capture the seasonal cycles of the vegetation mostly as a combination of the spectral signatures to derive in fine the

vegetation indices, which justifies the effort here.

S1 TOC (daily synthesis) products are inputs data valid for each pixel (available from http://www.vito-eodata.be). All land target are considered, snow and ice cover inclusive, at the exception of inner seas. The input data sets therefore resume to:

* TOC BLUE, RED, NIR at 1/3 km (VNIR) and SWIR at 2/3 km (SWIR)

* Radiometric quality flag for bands BLUE, RED, NIR, SWIR from Status Map (SM)

* Viewing and solar azimuth angles [0, 360°]

* Viewing Zenith Angles (VZA) and Viewing Azimuth Angles (VAA) are different for VNIR and SWIR

The latitudinal information is foremost ancillary information to estimate the value of the solar zenith angle. This latter corresponds here to noontime in the case of directional albedo and at 10 AM for TOC-R. The semi-empirical BRDF model of Roujean et al. (1992) further adjusted to the measurements reads:

$$R_j(\theta_{vj}, \theta_{sj}, \varphi_j) = \sum_{i=0}^{m-1} k_{lj} f_i(\theta_{vj}, \theta_{sj}, \varphi_j) \quad j = 1, n \quad (1)$$

where *n* represents the set of surface reflectance values, *m* stands for the number of kernels, equal to two here to account for the geometric and volume scattering kernels of Roujean et al. (1992), whereas $k_i$ are the model parameters and $f_i$ are the angular kernel functions. The angle $\phi$ is the relative azimuth between the directions of illumination and scanning. The angles $\theta_s$ and $\theta_v$ stand for the solar and viewing zenith angles, respectively. The algorithm first estimates the Black-Sky Albedo (BSA) or Directional-Hemispherical Reflectance (DHR) and White-Sky Albedo (WSA) or Bi-directional-Hemispherical Reflectance (BHR) in the four instrument channels by using the retrieved BRDF coefficients from Eq. (1). In the inversion procedure, a weight $w_j$ is assigned to each reflectance value, which is scaled as the inverse of the uncertainty on this reflectance value (e.g. Press et al., 1995). This uncertainty estimate is inherited from a statistical analysis of atmospherically corrected satellite scenes (see Geiger et al., 2008).

The solution to the linear least square inverse problem is stated as.

$$R = Fk. \quad (2)$$

The reflectance vector $b_j = R_j w_j$ is the solution of the following equation.

$$(A^T A)k = A^T b. \quad (3)$$

The design matrix $A_{ji} = F_{ji} w_j$ tailors the uncertainty covariance matrix:

$$C_k = (A^T A)^{-1}. \quad (4)$$

Provided few PROBA-V measurements and poor angular sampling scenario, more robust techniques like singular value decomposition (SVD) and QR-decomposition are advised to minimize numerical errors. The system can be better conditioned by adding constraints (e.g., Li et al., 2001; Hagolle et al., 2004; Pokrovsky et al., 2003). A priori information is given on the BRDF coefficients in terms of the first and second moments (average and standard deviation, respectively) of their a priori probability distribution function (PDF). That is:

$$k_i = k_{iap} \pm \sigma_{ap} [k_i] \quad (5)$$

Therefore, Eq. (3) is rewritten in the form.

$$(A^T A + c_{ap}^{-1})k = A^T b + c_{ap}^{-1} k_{ap} \quad (6)$$

The covariance matrix $C_{ap}$ for a priori information reduces to diagonal matrix terms by making the assumption of uncorrelated a priori information. Such matrix reads:

$$C_{ap} = C_k^{in} (1 + \Delta)^{(t_0 - t_{in}) / \Delta t} \quad (7)$$
<figure><img src="image_5.png" /><figcaption>Fig. 2. Effective temporal weight function in the recursive composition scheme.</figcaption></figure>

where the quantity $t_0 - t_{in}$ measures the elapsed time between two successive runs of the algorithm, with $\Delta t=10$ days. The characteristic time scale $\Delta$ is the full width at half mean of the weighing function (Fig. 2). Based on Eq. (7), the BRDF model inversion assigns less importance to observations as they get away from day $t_0$. Once the BRDF from the vector $k$ is known, the spectral albedo values are computerized by performing an angular integral of the kernel functions, $f_i$. Finally, a spectral integration of the spectral albedo leads to broadband albedos. For more about the algorithm, we will refer the reader to the implemented method for Spinning Enhanced Visible and InfraRed Imager (SEVIRI) onboard MSG (Meteosat Second Generation) as described in Geiger et al. (2008) from which the algorithm is basically derived.

## 3.2. Uncertainty estimate

The uncertainty estimate reports on the model performance by paying attention to the angular sampling and to the number of data. It is assumed that the probability distribution functions (PDF) of the errors of the TOC reflectance values are Gaussian and mutually uncorrelated. Using linear expressions, the albedo uncertainty estimate is obtained by propagating those of the model parameters as shown in Eq. (8) with the appropriate angular kernel integrals, $I$:

$$ \sigma[a] = \sqrt{I^T C_k I} \qquad (8) $$

Any residual cloud contamination will have an impact in the contribution of outliers in the PDF of the TOC reflectance errors, then on inversion results and surface albedo quality. To be outlined Eq. (8) remains valid for TOC-R.

## 3.3. Age of the product

The accumulation of information over days is mandatory for polar orbiting systems. The longer is the composite period, the more data to be stored for further exploitation. This is a trade-off process. Herein, PROBA-V land surface albedo and TOC-R products are reduced to 20-

<table><caption>Table 2
Albedo product quality flag information.</caption><thead><tr><th>Bit</th><th>Description</th><th>Value</th><th>Key</th></tr></thead><tbody><tr><td rowspan="4">Bits 0–1</td><td rowspan="4">Land Sea Mask</td><td>00</td><td>Sea</td></tr><tr><td>01</td><td>Land</td></tr><tr><td>10</td><td>Corrupted</td></tr><tr><td>11</td><td>Inland water</td></tr><tr><td>Bit 2</td><td>Input data</td><td>0</td><td>No</td></tr><tr><td rowspan="2">Bit 3</td><td rowspan="2">Majority Rule for snow or snow-free</td><td>1</td><td>Yes</td></tr><tr><td>0</td><td>No</td></tr><tr><td rowspan="2">Bit 4</td><td rowspan="2">A Priori Information</td><td>1</td><td>Yes</td></tr><tr><td>0</td><td>No</td></tr><tr><td rowspan="2">Bit 5</td><td rowspan="2">Snow</td><td>1</td><td>Yes</td></tr><tr><td>0</td><td>No</td></tr><tr><td rowspan="2">Bit 6</td><td rowspan="2">Aerosol correction</td><td>1</td><td>Yes</td></tr><tr><td>0</td><td>Climatology</td></tr><tr><td rowspan="3">Bit 7</td><td rowspan="3">Algorithm Failure</td><td>1</td><td>Near real time</td></tr><tr><td>0</td><td>No</td></tr><tr><td>1</td><td>Yes</td></tr></tbody></table>

days composite period against 30-days for 1 km VEGETATION program because the 333 m resolution of PROBA-V permits a better cloud removal. The product frequency is still 10 days. The synthetic period ends at days 10, 20, and last day of the month (Fig. 3). The true age of the product, Z_Age, is given to the user per pixel as a novelty. It is the median date value of clear sky dates over the composite period. Smoothness is slightly enhanced by assigning weight to the a priori information. Duration of cloudiness entails the propagation of a priori information and gap-filling.

## 3.4. Quality flag

The Quality Flag (Q-Flag) integrates all uncertainty assessments that led to the elaboration of the product (see Table 2). It reports merely information about a priori, snow coverage, and land pixel. If snow pixels are present during the composite period, the majority rule will apply. It means snow and snow-free pixels will compete in number to make an unmerged albedo product. The Bit 3 informs the user about the procedure application, with the risk for a merged albedo according to snow detection. The Bit 4 indicates that less than three observations were available. A priori estimate is reported in this case as a remedy to the gap-filling procedure. The Bit 6 is always 0 since a climatology is hitherto adopted for aerosol compound.

## 3.5. Narrow to broadband conversion

The narrow to broadband conversion is approximated as a weighted sum of the integrand at discrete values of the integration variable. The SWIR channel shows a coarser resolution compared to other channels and therefore a same SWIR pixel may be considered for neighboring pixels of other channels. Broadband albedo estimate at spectral interval $\gamma$ is derived from the spectral quantities by applying a linear transformation where $a_\gamma$ and $a_\beta$ stand for broadband and spectral albedos, respectively. Viz:

$$ a_{\gamma} = c_{0\gamma} + \sum_{\beta} c_{\beta\gamma} a_{\beta} \qquad (9) $$

<figure><img src="image_6.png" /><figcaption>Fig. 3. Illustration of the achievement of Z Age information.</figcaption></figure>
<table><caption>**Table 3**
Narrow- to broadband conversion coefficients for PROBA-V channels and uncertainty estimates.</caption><thead><tr><th>Broad-band</th><th>NDVI</th><th>c<sub>σγ</sub></th><th>c<sub>βγ</sub> Blue</th><th>c<sub>βγ</sub> Red</th><th>c<sub>βγ</sub> NIR</th><th>c<sub>βγ</sub> SWIR</th><th>RMSE</th></tr></thead><tbody><tr><td rowspan="2">Visible [0.4–0.7 µm]</td><td>[-1.,1.]</td><td>0.0010</td><td>0.5039</td><td>0.4923</td><td>–</td><td>–</td><td>0.0067</td></tr><tr><td>&lt; 0.2 &amp; snow</td><td>0.0284</td><td>0.5736</td><td>0.3837</td><td>–</td><td>–</td><td>0.0199</td></tr><tr><td rowspan="2">Near Infrared [0.7 - 4 µm]</td><td>[-1.,1.]</td><td>0.0140</td><td>–</td><td>0.0068</td><td>0.5677</td><td>0.3481</td><td>0.0135</td></tr><tr><td>&lt; 0.2 &amp; snow</td><td>0.0212</td><td>–</td><td>0.0438</td><td>0.5509</td><td>0.3633</td><td>0.0128</td></tr><tr><td rowspan="2">Total [0.3–4 µm]</td><td>[-1.,1.]</td><td>0.0097</td><td>0.1863</td><td>0.2212</td><td>0.3434</td><td>0.1817</td><td>0.0089</td></tr><tr><td>&lt; 0.2 &amp; snow</td><td>0.0248</td><td>0.1196</td><td>0.2764</td><td>0.3566</td><td>0.07</td><td>0.0154</td></tr></tbody></table>

<table><caption>**Table 4**
Uncertainty metrics for product validation.</caption><thead><tr><th>Gaussian Statistics</th><th>Comment</th></tr></thead><tbody><tr><td>Scatter plot of ground versus product</td><td>Qualitative assessment of agreement.</td></tr><tr><td>N: Number of samples</td><td>Indicative of the power of the validation</td></tr><tr><td>RMSE: Root Mean Square Error</td><td>RMSE computed between ground and product values should be compared to the RMSE value corresponding to ground measurements. Indicates the Accuracy (Total Error).</td></tr><tr><td>B: MeanBias</td><td>Relative values between the average of x and y were also computed. Difference between average values of ground and product. Indicative of accuracy and possible offset. Relative values between the average of x and y were also computed.</td></tr><tr><td>S: Standard deviation</td><td>Standard deviation of the pair differences. Indicates precision.</td></tr><tr><td>R<sup>2</sup>: Correlation coefficient</td><td>Indicates descriptive power of the linear accuracy test. Pearson coefficient was used.</td></tr><tr><td>Major Axis Regression (slope, offset)</td><td>Indicates some possible bias</td></tr><tr><td>p-Value</td><td>Test on whether the slope is significantly different to 1 (Null hypothesis: slope = 1). (p > 0.05 accepted)</td></tr><tr><td>% GCOS requirements</td><td>Percentage of pixels matching the GCOS requirements<sup>*</sup>.</td></tr></tbody></table>

<table><caption>**Table 5**
Summary of the Quality Assessment procedure.</caption><thead><tr><th>Quality Criteri</th><th>Product evaluated</th><th>Reference Product</th><th>Coverage</th></tr></thead><tbody><tr><td>Spatial Consistency</td><td>PROBA-V SA<br>Collection 300 m<br>Visual inspection of global maps.<br>Difference maps & global Scatter-plots (R<sup>2</sup>, RMSE, Bias, Scattering).</td><td>MODIS C6</td><td>European Region</td></tr><tr><td>Temporal Consistency</td><td>PROBA-V SA<br>Collection 300 m<br>Qualitative inspection of temporal variations.</td><td>MODIS C6 Ground measurements</td><td>EUVAL Validation sites</td></tr><tr><td>Overall Statistical Consistency</td><td>PROBA-V SA<br>Collection 300 m<br>Scatter-plots (R<sup>2</sup>, RMSE, Bias, Scattering) between pair of products.</td><td>MODIS C6</td><td>EUVAL (section 0)</td></tr><tr><td>Accuracy Assessment (Error)</td><td>PROBA-V SA<br>Collection 300 m<br>& MODIS C6<br>Scatter-plots, Pearson's correlation, Root Mean Square Error (RMSE), bias, linear fit (offset and slope, MAR).</td><td>Ground measurements</td><td>Validation sites</td></tr></tbody></table>

with coefficients $c_{σγ}$ and $c_{βγ}$ shown in Table 3. Three different broadband albedo intervals are considered: the total short-wave range from 300 nm to 4000 nm (BB), the visible wavelength range from 400 nm to 700 nm (VI), as well as the near infrared range from 700 nm to 4000 nm (NI). Negative NDVI values combining RED and NIR PROBA-V channels served to distinguish snow from snow-free pixels. The values in Table 3 were initially determined by van Leeuwen and Roujean (2002) for VEGETATION, further adapted to PROBA-V by VITO. Linear relationships were calibrated using synthetic data sets generated by the SAIL (Scattering by Arbitrarily Inclined Leaves) radiation transfer code (Verhoef, 1984), plus ASTER spectral library (Hook, 1998). They hold for both DHR (AL-DH) and BHR (AL-BH) products by assuming a spectral irradiance regardless to atmospheric effects. The errors of the spectral albedo estimates are supposed uncorrelated (see Table 3).

# 4. Tools of quality protocol for accuracy assessment

## 4.1. Information on product distribution

The processing chain is operated at the premises of VITO within the context of the EU Copernicus Global Land Service (CGLS). New land surface albedo and TOC-R products are proposed every ten days. They are distributed within three days after the last scene acquisition. For time being, they are available from the year 2014. The products are distributed as global files, combining the relevant layers into one single multi-band netCDF4-CF1.6 compliant format. To keep the file sizes to a reasonable size, not all layers are combined into one single file. Instead the TOC-R products are provided per band as well as the albedo products. The products can be downloaded from http://land.copernicus.eu/global and can be reformatted, cropped or subset through an on-line customization tool.

## 4.2. Methodology for quality assessment of the product

Preliminary Quality Assessment of the PROBA-V 1 Collection 300 m was focused on the surface albedo product merely. The evaluation is carried on over a European region covering the area from 35° to 65° of latitude and from -20° to 30° of longitude. The harnessing lumps together a cross-comparison with MODIS C6 products and then a direct validation using available ground observations coming from European Fluxes Data Cluster (EFDC) stations. The following main criteria were devised: spatial consistency, temporal consistency, and statistical assessment of discrepancies with similar products and accuracy. Surface albedo products receive regularly updated specifications (e.g. GCOS, 2016) as it is the case for an ECV (Essential Climate Variable). For surface albedo values lower than 0.05, the accuracy required is 0.0025. Therefore, a relative accuracy of 5% is still considered. Additional requirements come from the “WMO Rolling Requirement Review” that aids the setting of the priorities to be agreed by WMO (World Meteorological Organization) members and their space agencies. The GCOS requirements are only partly consistent with this process in that they provide only target but not “breakthrough” or “threshold” (i.e. minimum) requirements. But GCOS provides requirements on stability

<footnote>* The GCOS requirements on accuracy were used: Max(5%; 0.0025). Furthermore, an additional target level of Max(10%; 0.005) was used.</footnote>
<figure><img src="image_7.png" /><figcaption>Fig. 4. Location of the 109 selected European validation (EUVAL) sites and their classification into the main biome type. Red square delineates the region under study. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

<caption>**Table 6**</caption>

Selected 10 EFDC sites for Accuracy Assessment providing field albedo and diffuse measurements. The typical tower height used in this study is 10 m and then the footprint of the albedo measurements is 127 m (see Román et al., 2009).

<table><thead><tr><th>Code</th><th>Site name</th><th>Country</th><th>Land Cover</th><th>Lat (deg)</th><th>Lon (deg)</th><th>Diffuse Method</th></tr></thead><tbody><tr><td>CZ-BK1</td><td>Bily Kriz Forest</td><td>Czech Rep.</td><td>NLF</td><td>49.50472</td><td>18.54111</td><td>Direct</td></tr><tr><td>DE-Obe</td><td>Oberbärenburg</td><td>Germany</td><td>NLF</td><td>50.78362</td><td>13.71963</td><td>Indirect</td></tr><tr><td>ES-CPa</td><td>Cortes de Pallas</td><td>Spain</td><td>Shrublands</td><td>39.22417</td><td>-0.90305</td><td>Direct</td></tr><tr><td>ES-LMa</td><td>Majadas del Tietar</td><td>Spain</td><td>Savanna</td><td>39.9415</td><td>-5.77336</td><td>Direct</td></tr><tr><td>FR-Pue</td><td>Puechabon</td><td>France</td><td>EBF</td><td>43.74139</td><td>3.595833</td><td>Direct</td></tr><tr><td>IT-Col</td><td>Collelongo</td><td>Italy</td><td>BDF</td><td>41.84936</td><td>13.58814</td><td>Direct</td></tr><tr><td>IT-MBo</td><td>Monte Bondone</td><td>Italy</td><td>Herbaceous</td><td>46.01468</td><td>11.04583</td><td>Direct</td></tr><tr><td>IT-Tor</td><td>Torgnon</td><td>Italy</td><td>Shrublands</td><td>45.84444</td><td>7.578055</td><td>Indirect</td></tr><tr><td>DE-Tha</td><td>Tharandt</td><td>Germany</td><td>Mixed Forest</td><td>50.96361</td><td>13.56694</td><td>Direct</td></tr><tr><td>PL-Brd</td><td>Brody<sup>a</sup></td><td>Poland</td><td>Crop</td><td>52.43418</td><td>16.29952</td><td>Direct</td></tr></tbody></table>

<caption>**Table 7**</caption>

Additional 7 EFDC sites for temporal realism.

<table><thead><tr><th>Code</th><th>Site name</th><th>Country</th><th>Land Cover</th><th>Lat (deg)</th><th>Lon (deg)</th></tr></thead><tbody><tr><td>DE-Akm</td><td>Anklam</td><td>Germany</td><td>Shrublands</td><td>53.866170</td><td>13.683420</td></tr><tr><td>ES-ES</td><td>El Saler-Sueca</td><td>Spain</td><td>Crop</td><td>39.275550</td><td>-0.315278</td></tr><tr><td>IT-CA1</td><td>Castel d'Asso1</td><td>Italy</td><td>Crop</td><td>42.38041</td><td>12.02656</td></tr><tr><td>IT-CA2</td><td>Castel d'Asso2</td><td>Italy</td><td>Crop</td><td>42.37722</td><td>12.02604</td></tr><tr><td>IT-CA3</td><td>Castel d'Asso 3</td><td>Italy</td><td>Crop</td><td>42.38</td><td>12.0222</td></tr><tr><td>DE-Kli</td><td>Klingenberg</td><td>Germany</td><td>Crop</td><td>50.89288</td><td>13.52251</td></tr><tr><td>PL-Tuc</td><td>Tuczno</td><td>Poland</td><td>Wetland</td><td>53.192944</td><td>16.097472</td></tr></tbody></table>

that are not currently included in the WMO requirements database. The “WMO Observing Requirements Database” specifies requirements on the surface albedo for climatologic applications at three uncertainty quality levels: Goal (5%), Breakthrough (7%) and Threshold (10%).

Quality Assessment exercise follows here a procedure with pre- valent protocols and metrics defined to be consistent with the re- commendations of the Land Product Validation (LPV) group of the Committee on Earth Observation Satellite (CEOS) for the validation of satellite-derived land products. Spatial consistency mirrors the level of truth and repeatability of the spatial distribution of retrievals in the lack of spurious patterns or other artifacts (e.g., missing values, stripes, unrealistic low values, etc.). It is achieved through systematic visual

inspection analysis of the global maps, using the maps difference at a monthly basis. Temporal consistency reports on the degree of realism of the temporal variations. Temporal courses of satellite products were investigated over the confident sites of the EUVAL network and other selected stations. This supposes high frequency data are available like daily ground observations. Statistical analysis was performed over a selection of representative European sites for validation (EUVAL) based on uncertainties metrics associated to the scatter-plots between pairs of products (Table 4). To be outlined that only sites that are spatially re- presentative at the kilometer scale were put in valor. Table 5 sum- marizes the number of validation metrics used to verify the consistency between PROBA-V SA Collection 300 m for the whole 2014 year, using as references MODIS C6 product and ground data coming from EFDC stations.

A trimmed analysis of the spatial and statistical consistency was performed at 1 km resolution. PROBA-V SA 300 m products were re- sampled using an average value over a 3 × 3 pixels window (from 1/ 3 km to 1 km spatial resolution). MODIS C6 data were re-sampled using an average over 2 × 2 pixels window (from 1/2 km to 1 km). Twofold analysis based on temporal consistency and accuracy assessment relied on the native spatial resolution of each product. The comparison is carried on 10-days periods defined as the closest to the center of the temporal composite window, being 20 days (PROBA-V SA 300 m) and

<footnote><sup>a</sup> For PL-Brd Cropland site. Only Leaf-Off season (fall and winter) was considered for the accuracy assessment, due to the spatial homogeneity is not ensured during the Leaf-On season (spring and summer).</footnote>
<figure><img src="image_8.png" /><figcaption>Fig. 5. Maps of AL-DH (Left side) and AL-DH error estimates (right side) for several spectral domains and selected dates: AL-DH-BB for 2014.04.10 (Top), AL-DH-NI for 2014.06.10 (Middle) and AL-DH-VI for 2014.08.10 (Bottom).</figcaption></figure>

16 days (MODIS C6). Only pixels with Z_Age values beyond 10 days have been retained. This offers more realistic results of comparison by discarding products from sensors too far away in time. The date assigned to a PROBA-V pixel was established as the date given in the product file name - ending the composite window frame - less the Z_Age value.

## 4.3. Data sets for validation

### 4.3.1. PROBA-V data sets

The quality assessment is performed against a data set of PROBA-V SA V1 Collection 300 m images covering 15 tiles encompassing a wide European region from 35° to 65° of latitude and from -20° to 30° of longitude. The dataset is available at 10 days frequency covering the whole year 2014. For each tile, files in HDF5 format were considered with the content of black-sky and white-sky albedos in visible, NIR and

shortwave domains, plus ancillary information (error estimate, Quality Flag, Z_Age).

### 4.3.2. MODIS collection 6

For the cross-comparison, it is considered the MODIS BRDF/Albedo (MCD43A3) Collection 6 (DOI: https://doi.org/10.5067/MODIS/MCD43A3.006) at 500-meter spatial resolution, which achieved validation stage 3 in the CEOS/LPV hierarchy (Wang et al., 2018). Collection 6 provides improved quality compared with previous Collection 5, as well as more retrieval at high latitudes from use of all available observations. Only high quality MODIS products were retained for further analysis. It includes notably both directional hemispherical reflectance (black-sky albedo) at local solar noon and bi-hemispherical reflectance (white-sky albedo) for three broad-bands (visible: 0.3–0.7 µm, NIR: 0.7–5.0 µm, and Total: 0.3–5.0 µm). The MCD43A3 albedo quantities are disseminated on a sinusoidal grid, with temporal
<figure><img src="image_9.png" /><figcaption>Fig. 6. Maps of ALDH Z_Age for several dates. From Top to Bottom and from Left to Right: 2014.01.20, 2014.03.20, 2014.05.20, 2014.07.20, 2017.09.20 and 2014.11.20.</figcaption></figure>

composite period of 16 days. Both Terra and Aqua data are used in the generation of this product. As a quick remind, MODIS albedo algorithm uses atmospherically corrected reflectance data (MOD09 product with flags about “cloud”, “cirrus high” and “aerosol high”) to best fit the kernel-driven BRDF model Ross Thick Li Sparse-Reciprocal (RTLSR) (Lucht and Lewis, 2000). Only periods with 7 clear observations lead to full inversion. Otherwise it makes use of a backup algorithm with prior information for a magnitude inversion in the MODIS BRDF products. A gap filling occurs for less than three clear observations. To be outlined here that the MODIS albedo team provides separated snow-free gap filled products: https://www.umb.edu/spectralmass/terra_aqua_modis/v006/mcd43gf_cmg-gap filled_snow_free_products (Sun et al., 2017). Integrating BRDF model parameters lead to spectral albedos further converted into broadband albedos (Liang et al., 1999). A well-behaved agreement (r² = 0.82) was found between the comparison of mean yearly MODIS albedo Collection 5 retrievals with ground measurements

taken at 53 FLUXNET homogeneous sites (Cescatti et al., 2012).

### 4.3.3. Sites selection from European validation (EUVAL)

The temporal and statistical consistency was computed over 109 sites from the European Validation (EUVAL) network (Fig. 4). Only homogeneous sites covering at least 3 km² were deemed dependable based on previous study based on the biome type classification (GLC2000 (e.g. https://land.copernicus.eu/global/sites/cglsvito.be/files/products/GIOGL1_QAR_LAI300m-V1_11.10.pdf). The percentage of EUVAL sites per main biome is decomposed as follows: 14.37% for Broadleaf Forest, 21.56% for Needle-leaf forest, 32.68% for Cultivated areas, 9.8% for Shrublands, 11.76% for Herbaceous, and 9.80% for bare soils.

### 4.3.4. Accuracy assessment from EFDC

The accuracy assessment of PROBA-V SA Collection 300 m products
<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /></figure>

<figure><img src="image_13.png" /><figcaption>Fig. 7. AL-DH-VI difference maps (left side) and associated scatter plots (right side) between PROBA-V SA Collection 300 m and MODIS C6 products for the 20th March 2014 (Top) and 20th July 2014 (Bottom). Pixels with error estimates > 0.15 and Z_Age values > 10 days in cases of PROBA-V product were not considered for the analysis.</figcaption></figure>

finally prevails with stations from the European Fluxes Database Cluster (EFDC) (http://gaia.agraria.unitus.it/home/log-in/). These latter are coping with the mandatory criteria aforementioned. Land cover characteristics were obtained at the spatial resolution of 500 m based on Google Earth™ for sites matching first the requirement of homogeneity around the tower flux (Román et al., 2009). Both daily and noon measurements were harnessed. Different strategies have been applied to the ground measurements in order to fairly match with satellite data, which depends on the instrument. In case of PROBA-V, the value of the composite period is assigned to the true date. In case of MODIS C6, the composite period is 16 days and the date corresponds to the center of the 16-day composite period (9th day). A composite value is therefore built only if 70% at least of daily ground measurements are available. Table 6 displays the list of 10 homogeneous sites used in this study. Most of the EFDC stations measure the diffuse down-welling shortwave

radiation information. If not then estimate of diffuse radiation was performed in using the aerosol optical depth from the nearest AERONET station and the MODTRAN code. In addition, the ground information from 7 EFDC (Table 7) was only used to appraise the temporal shape as the diffuse fraction was not available.

# 5. Results of quality assessment

## 5.1. Product content

Spatially consistency of pattern distribution of surface albedo and associated error estimates are conspicuous from Fig. 5. The completeness of the mapping is a benefit from the recursive temporal scheme. Consistent values of Z_Age were found for the whole region during the year 2014, with values typically ranging between 5 and 15 (see Fig. 6).
<figure><img src="image_14.png" /></figure>

<figure><img src="image_15.png" /></figure>

<figure><img src="image_16.png" /><figcaption>Fig. 8. Idem Fig. 7 for AL-DH-NI.</figcaption></figure>

<figure><img src="image_17.png" /></figure>

As expected, higher Z_Age values were encountered for scenarios suffering from persistent cloud coverage (wintertime period and Northern latitudes).

## 5.2. Spatial consistency analysis

The set of Figs. 7–9 shows deviations in geographic distribution and scatter plots between PROBA-V Collection 300 m and MODIS C6 black-sky albedo products at mid-March and mid-July 2014. Snow flagged pixels were discarded from the analysis due to differences of strategy between the two sensors. For PROBA-V, rare pixels with error field > 0.15 were removed from the analysis as it is known possible residual cloud contamination with this first PROBA-V collection. Salient feature is the good spatial consistency that exists in visible domain (Fig. 7), with only differences of ± 0.025 in reflectance units. Also, dependable correlations (R² between 0.64 and 0.93) and mean bias close to zero are

evidenced. A slight tendency noticed with PROBA-V is lower values obtained for a low albedo (< 0.1), and higher values for a high albedo (> 0.1). This is a general trend observed if one excerpts the months of January and December 2014. As for NIR domain (Fig. 8), a good consistency is also reached with almost all differences within ± 0.05. The tendency is rather a random bias between PROBA-V SA Collection 300 m and MODIS C6, with opposite signs of the bias at some places for different periods. A good statistical rendering is obtained (R² around 0.8, RMSE between 0.02 and 0.03). Finally, the upshots of analysis for shortwave products (Fig. 9) closely follow the ones for NIR with nevertheless higher values compared to MODIS C6 products for all dates and regions again with the exception of northern latitudes and winter time situations.
<figure><img src="image_18.png" /></figure>

<figure><img src="image_19.png" /></figure>

<figure><img src="image_20.png" /></figure>

<figure><img src="image_21.png" /><figcaption>Fig. 9. Idem Fig. 7 for AL-DH-BB.</figcaption></figure>

5.3. *Temporal consistency analysis*

Temporal profiles of PROBA-V and MODIS C6 of black-sky and white-sky albedo products were analyzed for the 109 EUVAL sites and over 17 ground stations where ground measurements were available. Results are displayed for the three spectral ranges and the nominal spatial resolution of the product (333 m for PROBA V, 500 m for MODIS) (Fig. 10).

Vertical bars of PROBA-V Collection 300 m correspond to the un-
certainty assessment in the generation of the associated albedo product.
Blue-sky albedos are shown whenever daily ground observations are
available. Main lesson learned is that for Broadleaf Forests, PROBA-V
well captures the rapid seasonality of albedo from dense to low vege-
tation coverage, providing reliable and smooth trajectories. Only over
some sparse situations (a ratio of 25%), unexpected instabilities not
observed in reference products were observed in case of Broadleaf

Forests (see *Asturias* site during August–September). Although it is
prevalent in NIR, the shape also affects the shortwave band by spectral
construction. Over Needle-leaf Forests, unexpected variations (as
compared with reference products) in around 40% of cases were found
for some period (see B2.1#244 during July and *Cazorla NP* during
August–September). As for Cropland sites, close patterns are noticed
between sensors. To be outlined that PROBA-V well mimics the dy-
namic of surface albedo in link to phenological status. Such findings
were deemed consistent with MODIS C6 temporal variations and with
daily ground observations whenever available (see *El SalerSueca, Castel
d'Asso2*). Smoother transitions are conspicuous in the case of PROBA-V
due to the beneficial effect of a priori information that precludes dodgy
scenarios. Over Herbaceous and Shrublands, compliant temporal tra-
jectories are noticeable between sensor products and also ground daily
observations (e.g. *Anklam site* and *Semi-arid Almeria site*). Error esti-
mates were found generally concurring in magnitude with albedo
<figure><img src="image_22.png" /><figcaption>**Fig. 10.** Temporal profile of PROBA-V SA Collection 300 m (continuous line), and MODIS C6 (dashed line) AL-DH-BB (grey) and AL-DH-NI (yellow) and AL-DH-VI (blue) albedos over several selected sites. Vertical bars correspond to the associated error of PROBA-V SA Collection 300 m products. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

retrievals, excepted during the wintertime periods due to the occurrence of snow and grazing solar angles. MODIS C6 shows major gaps, which limits the comparison with PROBA-V. The snow episodes were correctly reported by PROBA-V in 70% of cases as compared to MODIS or ground daily observations (e.g. B2.1#259). Unreported snow events may explain by the choice of the majority rule to segregate snow or

snow-free pixels over the spanning composite window.

## 5.4. Blue-sky albedo

Blue-sky albedo is a genuine product from the point of view of in situ measurement, although not yet included in the Copernicus GLS.
<figure><img src="image_23.png" /></figure>

<figure><img src="image_24.png" /></figure>

<figure><img src="image_25.png" /></figure>

<figure><img src="image_26.png" /></figure>

<figure><img src="image_27.png" /></figure>

(caption on next page)
<figure><img src="image_28.png" /><figcaption>Fig. 12. Black-sky (Top) and White-Sky (Bottom) PROBA-V SA Collection 300 m versus MODIS C6 albedo products scatter-plots over all EUVAL sites during the 2014 year for PAR (left), NIR (center) and shortwave (right) domains. The terms B and S represent the mean and the standard deviation of the difference. Continuous black line corresponds to 1:1 line and dashed lines to GCOS uncertainty levels. Red line corresponds to the Major Axis Regression (MAR). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

Over stationary targets like forests (e.g. Billy Kriz Forest, Tharandt), PROBA-V products seem efficient for the reproduction, ditto MODIS C6 (Fig. 11). PROBA-V albedo shows however a slight seasonality in the bias, with larger values during the summer season. Over sites marked by topography (Monte Bondone, Torgnon), PROBA-V better agrees with

ground data. The surges in albedo values due to snow events in Monte Bondone and Torgnon were well captured by PROBA-V SA product. The case of Brody sites with crops also reveals the performance of PROBA-V product (Fig. 11).
<caption>**Table 8**</caption>

Relevant statistics between PROBA-V SA Collection 300 m versus MODIS C6 Albedos products over EUVAL sites for the 2014 year. *p*-value corresponds to the test on whether the slope is significantly different to 1.

<table><caption>PROBA-V SA Collection 300 m vs MODIS C6</caption><thead><tr><th></th><th>AL-DH-VI</th><th>AL-DH-NI</th><th>AL-DH-BB</th><th>AL-BH-VI</th><th>AL-BH-NI</th><th>AL-BH-BB</th></tr></thead><tbody><tr><td>Correlation (R<sup>2</sup>)</td><td>0.5</td><td>0.6</td><td>0.58</td><td>0.54</td><td>0.72</td><td>0.71</td></tr><tr><td>Bias</td><td>-0.007</td><td>-0.006</td><td>0.005</td><td>0.006</td><td>-0.002</td><td>0.013</td></tr><tr><td>RMSE</td><td>0.033</td><td>0.045</td><td>0.033</td><td>0.031</td><td>0.038</td><td>0.03</td></tr><tr><td>Offset (MAR)</td><td>0</td><td>-0.0</td><td>-0.02</td><td>0.01</td><td>0.001</td><td>0.01</td></tr><tr><td>Slope (MAR)</td><td>0.95</td><td>1.07</td><td>1.14</td><td>0.92</td><td>0.95</td><td>1.04</td></tr><tr><td>p-Value</td><td>&lt; 0.001</td><td>&lt; 0.001</td><td>&lt; 0.001</td><td>&lt; 0.001</td><td>&lt; 0.001</td><td>&lt; 0.001</td></tr><tr><td>%optimal (GCOS)</td><td>12</td><td>5.3</td><td>6.6</td><td>9.7</td><td>6.3</td><td>6.8</td></tr></tbody></table>

<figure><img src="image_29.png" /></figure>

<figure><img src="image_30.png" /><figcaption>**Fig. 13.** Accuracy Assessment of PROBA-V SA Collection 300 m (left) and MODIS MCD43A3 C6 (right) blue-sky albedo satellite products versus ground values coming from European stations during the 2014 year for Snow Free conditions. The terms B and S represent the mean and the standard deviation of the difference. Continuous black line corresponds to 1:1 line and dashed lines to GCOS uncertainty levels. Red line corresponds to the Major Axis Regression (MAR). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

### 5.5. Overall statistical consistency and accuracy assessment

The overall consistency of the PROBA-V SA Collection 300 m with reference MODIS C6 is resumed here with EUVAL sites for year 2014 (Fig. 12, Table 8).

For visible range, good overall consistency between PROBA-V and

MODIS C6 was found, with RMSE values around 0.03. It means PROBA-V is less than MODIS for low albedo values, mainly black-sky albedo. Slight negative bias (−0.007, 11%) was found for black-sky albedo, whereas the opposite trend (positive bias of 0.006, 9%) was observed for white-sky albedo. For the NIR, very good consistency was found, with RMSE values around 0.04 and very low negative mean bias (2.4% for AL-DH and 0.8% for AL-BH). Similar good results were found for Broad Band, showing low RMSE values of around 0.03 and low positive mean bias (< 8%). Different slopes were observed. This is reliant to the spectral albedo: slope < 1 for PAR, slope > 1 for Broad Band, with different trend in the slope for NIR. Opposite trends were observed depending on the albedo type (slope > 1 for AL-DH and slope < 1 for AL-BH). In all cases, slopes significantly differed to 1 with statistical significance (*p*-value < 0.05). To go further in the investigation of the appropriateness of the satellite albedo products (diagnosing “blue-sky” albedo), scatter plots versus field measurements were produced for the whole 2014 year over 10 EFDC sites of different vegetation types (Fig. 13, Table 9). Note that the number of samples (N = 286 in case of PROBA-V against N = 228 in case of MODIS C6) indicates the better completeness of PROBA-V, as already observed on the temporal profiles (Fig. 11), on which PROBA-V provides large number of retrievals. Overall accuracy of RMSE = 0.032 was found for PROBA-V SA Collection 300 m products, showing the same result in terms of RMSE as MODIS C6. Slight positive mean bias of 0.01 (9%) of PROBA-V SA Collection 300 m was found, showing improved results in case of MODIS C6, with no bias and relationship very close to the 1:1 line. In spite of positive results, the low percentage of pixels within GCOS requirements (5.9% in case of PROBA-V and 9.2% in case of MODIS C6) outlines the difficulty to achieve the requirements.

### 5.6. TOC-R product

Alongside developments on surface albedo, the normalized reflectance TOC-R appears as a surrogate for any application not necessarily reliant on the radiation transfer and energy budget. Such product targets a broad category of users eager to consider PROBA-V radiometry to catch up the phenology. A product like TOC-R is difficult to validate because it relies to a specific geometry serving for reference that should be obtained from ground using narrow field of view. In the case of PROBA-V, the TOC-R product required by users corresponds to a reflectance that would be measured at 10 AM local time. In practice, such virtual reflectance is obtained by using the adjusted BRDF coefficients to simulate it. In the case of MODIS, the strategy is different since the distributed normalized reflectance is a nadir BRDF-Adjusted Reflectance (NBAR). Furthermore, PROBA-V and MODIS spectra are not matching. Therefore, the option herein was to rather cross-compare NDVI products elaborated from TOC-R (PROBA-V) or NBAR (MODIS). To be noticed that the spatial resolution of PROBA-V is 300 m against 500 m for MODIS. Fig. 14 shows for the whole year 2014 the NDVI temporal profiles for varied land units belonging to different locations (France, Spain, Italy, UK, Tunisia).
<caption>**Table 9**</caption>

Relevant statistics of the Accuracy Assessment of PROBA-V SA Collection 300 m (left) and MODIS MCD43A3 C6 (right) blue-sky albedo satellite products versus ground values coming from European stations during the 2014 year for Snow Free conditions. p-value corresponds to the test on whether the slope is significantly different to 1.

<table><thead><tr><th rowspan="2"></th><th colspan="2">EFDC SITES (2014)</th></tr><tr><th>PROBA-V SA Collection 300 m</th><th>MCD43A3 C6</th></tr></thead><tbody><tr><td>Correlation (R<sup>2</sup>)</td><td>0.55</td><td>0.60</td></tr><tr><td>Bias</td><td>0.011</td><td>0.000</td></tr><tr><td>RMSE</td><td>0.032</td><td>0.032</td></tr><tr><td>Offset (MAR)</td><td>0.01</td><td>0.00</td></tr><tr><td>Slope (MAR)</td><td>1.01</td><td>0.99</td></tr><tr><td>p-value</td><td>0.784</td><td>0.919</td></tr><tr><td>%optimal (GCOS)</td><td>5.9</td><td>9.2</td></tr></tbody></table>

Twofold PROBA-V NDVI profiles are shown whenever considering the end of the composite period given in file name or accounting for Z_Age information. Significant shift in time is noticed. The cross-comparison with MODIS remains qualitative as long it is not the same spatial resolution. The cultivated site of Barrax (Spain) mirrors the same seasonal patterns for all plotted NDVI, with closest profile to MODIS C6 using PROBA-V Z_Age. Missing data are conspicuous for MODIS due to discarding data based on Quality Assurance. The NDVI of Chimbolton (UK) is garbled by cloudiness. But despite two sensors (TERRA and AQUA) against one sensor for PROBA-V, MODIS displays less and more erratic NDVI products. The site of Merguellil (Tunisia) with crops reveals a good correlation between MODIS and PROBA-V (using Z_Age) except at the onset where a bias is noticed. MODIS NDVI is above PROBA-V during winter months. Same conclusion applies for the deciduous forest of Harth (France) and the rice crop of Albufera (Spain). This remains valid for the cropped area of Castellaro (Italy) where PROBA-V NDVI shows at the onset a first vegetation peak not clearly depicted by MODIS C6. It follows shared features for a plateau and the decay.

## 6. Summary and future prospects

An initiative was undertaken to provide operational global estimates of the surface albedo (SA) and TOC-R 300 m and 1 km fields in the framework of the Global Land Service (GLS) of Copernicus. It is represented by the PROBA-V Collection 1. From the outcomes of the present study, it clearly offers a great potential in agricultural area for instance. It is glimpsed that product would improve owing to a better cloud detection, particularly ice cloud removal - at the contrary of MODIS no thermal band exists for PROBA-V - and aerosol signal correction. The occurrence to have ice cloud flagged as snow was noticed and can at least be solved by a temporal filter. Besides, the analysis was performed for snow free conditions because of the different strategy of temporal sampling between the two sensor systems. Therefore, the cross-comparison was limited to data sets ingestion beyond a value of 0.5 typically. The narrow to broadband conversion process also represents a possible source of future improvement for PROBA-V.

Time and space good consistency are noticed between the two sensors, with overall discrepancies (RMSE) of 0.03 for PAR and short-wave domains, and around 0.04 for NIR. PROBA-V Collection 300 m displayed reliable temporal trajectories, matching with temporal trend from MODIS C6 and daily ground observations. Smooth transitions may prevail as well as surges in albedo values, to be reasonably reproduced. The accuracy assessment performed over 10 EFDC homogeneous targets shows overall uncertainty (RMSE) of 0.032 with a mean bias of 0.01, which is about the same order as MODIS C6. Although this new production presents some spurious variability on short time scales, the comparison against MODIS C6 is deemed quite promising at this stage. It will even offer enhanced perspectives with next PROBA-V

<figure><img src="image_31.png" /></figure>

<figure><img src="image_32.png" /></figure>

<figure><img src="image_33.png" /></figure>

<figure><img src="image_34.png" /></figure>

<figure><img src="image_35.png" /></figure>

<figure><img src="image_36.png" /><figcaption>**Fig. 14.** Temporal profiles of PROBA-V TOC-R Collection 300 m using the date of the end of the composite period given in file name (blue) and the date corrected from Z_Age (red). The Nadir BRDF-Adjusted Reflectance MODIS C6 (MCD43A4) is reported (green). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

Collections, owing to users feedbacks. The TOCR-R product was analyzed through NDVI for cross-comparing with MODIS. The Z_Age field improves the matching between the two sensors. Less number of outliers is noticed for PROBA-V due to the recurrent strategy applied. The existing bias with MODIS for the wintertime quiescent period is due to the use of a slant path as reference in the case of PROBA-V, which seems
to offer more dynamic on the phenology.

## Acknowledgements

"This work was initiated under the FP7 Imagines project and thereafter continued and delivered through the land service of Copernicus, the Earth Observation program of the European Commission. The research leading to the current version of the product has received funding from various European Commission Research and Technical Development programs. The product is based on PROBA-V 300m data ((c) ESA and distributed by VITO)."

## References

Baret, F., Hagolle, O., Geiger, B., Bicheron, P., Miras, B., Huc, M., Berthelot, B., Weiss, M., Samain, O., Roujean, J.-L., Leroy, M., 2007. LAI, FAPAR and fCover CYCLOPES global products derived from VEGETATION. Part 1: principles of the algorithm. Remote Sens. Environ. 110, 275-286.

Barnsley, M.J., Strahler, A.H., Morris, K.P., Muller, J.P., 1994. Sampling the surface bidirectional reflectance distribution function (BRDF): evaluation of current and future satellite sensors. Remote Sens. Rev. 8, 271-311.

Cescatti, A., Marcolla, B., Vannan, S.K.S., Pan, J.Y., Roman, M.O., Yang, X., Ciais, P., Cook, R.B., Law, B.E., Matteucci, G., Migliavacca, M., Moors, E., Richardson, A.D., Seufert, G., Schaaf, C.B., 2012. Intercomparison of MODIS albedo retrievals and in situ measurements across the global FLUXNET network. Remote Sens. Environ. 121, 323-334.

Dierckx, W., Sterckx, S., Benhadj, I., Livens, S., Duhoux, G., Van Achteren, T., ... Saint, G., 2014. PROBA-V mission for global vegetation monitoring: standard products and image quality. Int. J. Remote Sens. 35 (7), 2589-2614.

GCOS-154, 2016. Systematic Observation Requirements for Satellite-based Products for Climate Supplemental Details to the Satellite-based Component of the Implementation Plan for the Global Observing System for Climate in Support of the UNFCCC - 2016 Update. WMO, Geneva, Switzerland.

Geiger, B., Carrer, D., Franchistéguy, L., Roujean, J.L., Meurey, C., 2008. Land surface albedo derived on a daily basis from Meteosat Second Generation observations. IEEE Transactions on Geoscience and Remote Sensing 46 (11), 3841-3856.

Hagolle, O., Lobo, A., Maisongrande, P., Cabot, F., Duchemin, B., de Pereyra, A., 2004. Quality assessment and improvement of temporally composited products of remotely sensed imagery by combination of VEGETATION 1 & 2 images. Remote Sens. Environ. 94, 172-186.

Hook, S.J., 1998. ASTER Spectral Library. https://speclib.jpl.nasa.gov/.

Hu, B., Lucht, W., Li, X., Strahler, A.H., 1997. Validation of kernel-driven models for global modeling of bidirectional reflectance. Remote Sens. Environ. 62, 201-214.

Jin, Y., Gao, F., Shaaf, C.B., Li, X., Strahler, A.H., Bruegge, C.J., Martonchik, J.V., 2002. Improving MODIS surface BRDF/albedo retrieval with MISR multiangle observations. IEEE Trans. Geosci. Remote Sens. 40, 1593-1604.

Justice, C.O., et al., 1998. The moderate resolution imaging spectroradiometer (MODIS): land remote sensing for global change research. IEEE Trans. Geosci. Remote Sens. 36, 1228-1249.

van Leeuwen, W., Roujean, J.-L., 2002. Land surface albedo from the synergistic use of polar (EPS) and geo-stationary (MSG) observing systems. An assessment of physical uncertainties. Remote Sensing of Environment 81, 273-289.

Leroy, M., Deuzé, J.L., Bréon, F.M., Hautecoeur, O., Herman, M., Buriez, J.C., Tanré, D., Bouffiès, S., Chazette, P., Roujean, J.L., 1997. Retrieval of atmospheric properties and surface bidirectional reflectances over the land from POLDER/ADEOS. J. Geophys. Res. 102 (D14), 17023-17037.

Li, X., Gao, F., Wang, J., Strahler, A., 2001. A priori knowledge accumulation and its application to linear BRDF model inversion. J. Geophys. Res. 106 (D11), 11925-11935.

Liang, S., Strahler, A.H., Walthall, C., 1999. Retrieval of land surface albedo from satellite observations: a simulation study. J. Appl. Meteorol. 38, 712-725.

Lucht, W., Lewis, P., 2000. Theoretical noise sensitivity of BRDF and albedo retrieval from the EOS-MODIS and MISR sensors with respect to angular sampling. Int. J. Remote Sens. 21 (1), 81-98.

Lucht, W., Roujean, J.L., 2000. Considerations in the parametric modeling of BRDF and albedo from multiangular satellite sensor observations. Remote Sens. Rev. 18, 343-379.

Muller, J.-P., et al., 2012. The ESA GlobAlbedo Project for Mapping the Earth's Land Surface Albedo for 15 Years From European Sensors. Paper Presented at IEEE Geoscience and Remote Sensing Symposium (IGARSS) 2012. IEEE, Munich, Germany (22-27.7.12).

Pokrovsky, I.O., Pokrovsky, O.M., Roujean, J.-L., 2003. Development of an operational procedure to estimate surface albedo from the SEVIRI/MSG observing system in using POLDER BRDF measurements. Remote Sens. Environ. 87, 198-242.

Press, W.H., Teukolsky, S.A., Vetterling, W.T., Flannery, B.P., 1995. Numerical Recipes in Fortran. Cambridge University Press.

Rahman, H., Dedieu, G., 1994. SMAC: a simplified method for the atmospheric correction of satellite measurements in the solar spectrum. Int. J. Remote Sens. 15, 123-143.

Román, M.O., Schaaf, C., Woodcock, C., et al., 2009. The MODIS (collection V005) BRDF/albedo product: assessment of spatial representativeness over forested landscapes. Remote Sens. Environ. 113 (11).

Roujean, J.-L., Leroy, M., Deschamps, P.-Y., 1992. A bidirectional reflectance model of the Earth's surface for the correction of remote sensing data. J. Geophys. Res. 97 (D18), 20455-20468.

Samain, O., Geiger, B., Roujean, J.-L., 2006. Spectral normalization and fusion of optical sensors for the retrieval of BRDF and albedo: application to VEGETATION, MODIS and MERIS data sets. IEEE Trans. Geosci. Remote Sens. 44 (11, Part 1), 3166-3179.

Schaepman-Strub, G., Schaepman, M.E., Painter, T.H., Dangel, S., Martonchik, J.V., 2006. Reflectance quantities in optical remote sensing—definitions and case studies. Remote sensing of environment 103 (1), 27-42.

Sterckx, S., Benhadj, I., Duhoux, G., Livens, S., Dierckx, W., Goor, E., ... Nackaerts, K., 2014. The PROBA-V mission: image processing and calibration. International journal of remote sensing 35 (7), 2565-2588.

Strahler, A.H., 1994. Vegetation Canopy Reflectance Modeling - Recent Developments and Remote Sensing Perspectives. Proceedings of the 6th International Symposium on Physical Measurements and Signatures in Remote Sensingpp. 593-600.

Strahler, A.H., Muller, J.P., et al., 1999. MODIS BRDF/Albedo Product: Algorithm Theoretical Basis Document, Version 5.0.

Stroeve, J.C., Nolin, A.W., 2002. New methods to infer snow albedo from the MISR instrument with applications to the Greenland ice sheet. IEEE Trans. Geosci. Remote Sens. 40, 1616-1625.

Sun, Q., Wang, Z., Li, Z., Erb, A., Schaaf, C.B., 2017. Evaluation of the global MODIS 30 arc-second spatially and temporally complete snow-free land surface albedo and reflectance anisotropy dataset. Int. J. Appl. Earth Obs. Geoinf. 58, 36-49.

Verhoef, W., 1984. Light scattering by leaf layers with application to canopy reflectance modeling, the SAIL model. Remote Sens. Environ. 16, 125-141.

Wang et al., 2018. Capturing rapid land surface dynamics with collection V006 MODIS BRDF/NBAR/albedo (MCD43) products. Remote Sens. Environ. 207, 50-64.

Wanner, W., Li, X., Strahler, A.H., 1995. On the derivation of kernels for kernel-driven models of bidirectional reflectance. J. Geophys. Res. 100 (D10), 21077-21090.

Wanner, W., Strahler, A., Hu, B., Lewis, P., Muller, J.-P., Li, X., Barker-Schaaf, C., Barnsley, M., 1997. Global retrieval of BRDF and albedo over land from EOS MODIS and MISR data: theory and algorithm. J. Geophys. Res. 102 (D14), 17143-17161.

Wolters, E., Dierckx, W., Iordache, M.-D., Swinnen, E., February 2017. PROBA-V products user manual, version 2.1. http://habistat.vgt.vito.be/sites/default/files/PROBAV-Products_User_Manual_v2.1.pdf.