Contributors: G. Schwaizer (ENVEO IT GmbH), T. Nagler (ENVEO IT GmbH), M. Heinrich (ENVEO IT GmbH), P. Malcher (ENVEO IT GmbH)

Issued by: ENVEO IT GmbH / Gabriele Schwaizer

Date: 30 Apr 2026

Ref: C3S2_WP3-DDP-SCE-01_202506_ATBD_v1.0

Official reference number service contract: 2024/C3S2_313d_ENVEO/SC1

<style>
.special_indent ul > li > ul {
    padding-left: 0;
}

.special_indent ul > li > ul > li > ul {
    padding-left: 1.5em;
}
</style>

History of modifications

List of datasets covered by this document

Acronyms

General definitions

Climate Data Record (CDR): Defines a time series of measurements of sufficient length, consistency, and continuity to determine climate variability and change.

Essential Climate Variable (ECV): a physical, chemical or biological variable or a group of linked variables that critically contributes to the characterization of Earth's climate.

Interim Climate Data Record (ICDR): Defines a dataset that has been forward processed, using the baselined Climate Data Record algorithm and processing environment but whose consistency and continuity have not been verified. Eventually, it will be necessary to perform a new reprocessing of the CDR and ICDR parts together to guarantee consistency, and the new reprocessed data record will replace the old CDR.

Snow Cover Extent (SCE): Areal extent of snow-covered land, which can be expressed as binary or as a fraction.

Snow Cover Fraction (SCF): Fraction of snow covered area per pixel, given in per cent.

L0: Unprocessed instrument and payload data at full resolution.

L1A: Reconstructed unprocessed instrument data at full resolution, time referenced, and annotated with ancillary information, including radiometric and geometric calibration coefficients and georeferencing parameters, computed and appended, but not applied, to L0 data.

L1B: Level 1A data that have been processed to sensor units.

L2: Retrieved environmental variables at the same resolution and location as the level 1 source.

L3C : L3 Collated: Level 2 variables mapped on a defined grid with reduced requirements for ancillary data. Observations combined from a single instrument into a space-time grid. A typical SCE L3C product may contain all the observations from a single instrument in a 24-hour period.

Executive summary

This document is the Algorithm Theoretical Basis Document (ATBD) for the Snow Cover Extent (SCE) Climate Data Record (CDR) and Interim CDR (ICDR), version 1.0. The SCE products are produced by ENVEO-Environmental Earth Observation IT GmbH using satellite data from NOAA-7/9/11/12/14 AVHRR/2, Terra MODIS, and Sentinel-3 A&B SLSTR sensors. The Terra MODIS-based SCE products are generated within the ESA Climate Change Initiative Extension (CCI+) Snow project, while the NOAA AVHRR/2 and Sentinel-3 SLSTR-based SCE products are produced under the Copernicus Climate Change Service (C3S). 

The ATBD describes the physical and mathematical basis of algorithms and systems used to generate the SCE products. It provides an overview of the satellite missions and instruments used for the generation of the SCE products, each contributing to a daily global SCE time series from 1982 to the present, with optimal spatial resolution.

From 1982 to 2000, data from NOAA missions allow SCE retrieval with a pixel spacing of 0.05 deg x 0.05 deg. Since late February 2000, the Terra MODIS sensor has enabled daily global SCE retrievals at an enhanced spatial resolution of 0.01 deg x 0.01 deg. The Sentinel-3 SLSTR satellite system, used from 2023 onwards, ensures the continuation of these spatial characteristics into the future.

The ATBD describes the satellite input and auxiliary datasets used to generate daily global SCE products. The retrieval process is modular, beginning with the pre-processing of satellite data, followed by cloud masking. For cloud-free pixels, a pre-classification is applied to identify areas likely to be snow-free. For all remaining pixels, the snow cover fraction per pixel is estimated. In forested areas, snow cover fraction viewable from above (SCFV) is distinguished from the snow cover fraction on the ground (SCFG), retrieved by applying a canopy correction. In non-forested regions, SCFV and SCFG are the same. 

Sensor-specific post-processing steps are applied. For AVHRR/2, the classification scheme is complemented by additional criteria tailored to the sensor's characteristics. For MODIS and SLSTR, tile-based classifications are merged into global daily snow cover fraction products. Daily integrated temperature and precipitation data are used to refine the cloud mask across all sensors. 

The ATBD provides a detailed description of the snow classification dependencies, the use and source of satellite and auxiliary data, the selected SCE retrieval algorithm and the associated uncertainty estimation, and a summary of the generated output products.

Missions and Instruments 

This section provides an overview on the missions and the relevant instruments used for the generation of the ECV Snow Cover Extent products. An overview on the temporal coverage per mission and the usage per mission for the generation of the ECV Snow Cover Extent products is provided in Table 1.1 and is illustrated in the Product User Guide and Specification document.

Table 1.1: Missions and relevant instruments used for the generation of ECV Snow Cover Extent products.

NOAA AVHRR

The Polar-orbiting Operational Environmental Satellite (POES) constellation of weather satellites is a joint effort between the National Oceanic and Atmospheric Administration (NOAA) and the European Organisation for the Exploitation of Meteorological Satellites (EUMETSAT). A series of platforms, with the first one launched in 1978, operate in both morning and afternoon sun synchronous low-Earth polar orbits (804 - 870 km, depending on the satellite). Each platform carries an Advanced Very High Resolution Radiometer (AVHRR) on board. AVHRR is a cross-track scanning system.

Information about the POES Radiometer data can be found at https://www.ncei.noaa.gov/access/metadata/landing-page/bin/iso?id=gov.noaa.ncdc:C00683 and https://www.eumetsat.int/our-satellites/metop-series.

Details about the NOAA AVHRR data used for the Snow Cover Extent retrieval are provided in Section 2.

Terra MODIS

Terra is an international mission launched by NASA in 1999, flying in low Earth orbit (705 km altitude). It carries five scientific instruments, one of which is the Moderate Resolution Imaging Spectroradiometer (MODIS). Terra MODIS data are used to generate the daily global Snow Cover Extent products from 2000 to 2022. 

Information about all instruments on board of Terra can be found at https://terra.nasa.gov/about/terra-instruments.

Details about the Terra MODIS data used for the Snow Cover Extent retrieval are provided in Section 2.

Sentinel-3 SLSTR

Sentinel-3 is an operational mission that is part of the Copernicus programme. There are currently two satellites in orbit, Sentinel-3A and Sentinel-3B, flying in constellation in low Earth orbit (800 - 830 km altitude). Sentinel-3A was launched in 2016, Sentinel-3B was launched in 2018. Each satellite carries four instruments on board, including the Sea and Land Surface Temperature Radiometer (SLSTR), a dual-view scanning temperature radiometer used for the extension of the Snow Cover Extent Climate Data Record from 2023 onwards.

Information about each instrument can be found at https://www.esa.int/Applications/Observing_the_Earth/Copernicus/Sentinel-3/Instruments.

Details about the Sentinel-3 SLSTR data used for the Snow Cover Extent retrieval are provided in Section 2.

Input and auxiliary data

This section provides detailed information on the satellite input data used for the generation of daily global SCE products. These include:

• NOAA AVHRR/2
• Terra MODIS
• Sentinel-3A&B SLSTR

Spectral bands used as input for the Snow Cover Extent product generation are listed in Table 2.1. Further details on the satellite data are provided in Section 2.1.

Table 2.1: Spectral bands per instrument used for the generation of ECV Snow Cover Extent products. 

*⁾ Reflective part is extracted as described e.g. Baum and Trepte (1999) using the parameterisation of Trishchenko (2006).

**) Band available but not used.

***) Band used only for cloud screening from MODIS and SLSTR data.

In addition, all auxiliary datasets used in the generation of daily global SCE products are described in detail. These include:

• C3S land cover classification map of the year 2000, used to identify different surface categories, support the preparation of other auxiliary layers, and mask water, and permanent snow and ice areas in the final products.
• Digital elevation model used as input for the Normalized Difference Snow Index (NDSI) threshold map (see next bullet). 
• Threshold map for the Normalized Difference Snow Index (NDSI), used to pre-classify areas being likely snow free, 
• Two-way forest transmissivity map, used for the canopy correction to estimate in forested areas the snow cover fraction on ground. 
• Sensor-specific spectral reflectance maps used for the snow cover fraction classification per pixel.

These auxiliary data sets are considered as constant over time, and do not reflect any changes in surface conditions.

Satellite data

AVHRR NOAA (EUMETSAT)

MODIS Terra (NASA)

SLSTR Sentinel-3 (Copernicus)

Auxiliary data

Land Cover Classification

The land surface classification for the year 2000, available from the C3S Land Cover (LC) classification (https://cds.climate.copernicus.eu/datasets/satellite-land-cover?tab=overview), is utilized as an auxiliary input within the processing chain. The year 2000 was selected as it approximately represents the midpoint of the temporal coverage of the satellite sensors employed in generating the Snow Cover Extent (SCE) Climate Data Record (CDR). This dataset is hereafter referred to as LC2000. 

 Water and permanent snow and ice areas

The LC2000 dataset is employed to mask water bodies and permanent snow and ice areas in the final SCE products. Specifically, water bodies (LC2000 class code 210) and permanent snow and ice areas (LC2000 class code 220) are extracted as separate binary masks (values 0 and 1), preserving the original spatial resolution of approximately 300 meters and the geographic (latitude/longitude) projection. These binary layers are subsequently aggregated to match the spatial resolutions of 0.05 deg. and 0.01 deg., corresponding to the SCE products derived from the AVHRR, MODIS, and SLSTR sensors, respectively. Aggregation is performed by computing the mean value of all mask pixels within each output grid cell, resulting in fractional layers that represent the proportion of water or permanent snow and ice area per pixel at the target resolution. Based on the specific pixel spacing of the sensor used for SCE product generation, thresholding is applied to these fractional layers to mask relevant pixels as water or permanent snow and ice.

 Salt lakes

In the LC2000 dataset, salt lakes are neither represented as a distinct land cover class nor included within the "water bodies" category. However, many salt lakes may become partially or fully inundated during periods of precipitation. Since water exhibits spectral characteristics similar to snow, applying a comprehensive water mask helps to minimize commission errors in the snow cover maps. Additionally, the dry surfaces of salt lakes—typically covered with a salt crust—can produce high values in the Normalized Difference Snow Index (NDSI) and high reflectance values in the visible spectrum, which are both key inputs in the SCE retrieval algorithm.

As no evidence of actual snow cover on salt lakes was identified, these areas were masked as a separate class to prevent misclassification. To this end, salt lakes were manually delineated in regions where misclassification was observed in the ESA CCI Snow Cover Fraction on Ground (SCFG) CRDP v2.0 product, primarily due to water or salt-related spectral confusion. Mapping focused on regions with a high occurrence of salt lakes, including Utah (USA), the Central and Southern Andes, Central Iran, Western Uzbekistan, and Southern Australia. More than 100 salt lakes were mapped, with the largest examples being Lake Eyre (Australia), Salar de Uyuni (Bolivia), and the Great Salt Lake (USA), each covering areas of approximately 10,000 km².

Forest classes

Selected classes of the LC2000 data set, listed in Table 2.2, are used for the preparation of the transmissivity map (cf. Section 2.2.4) to account for the shading effect of the forest canopy in the estimation of the snow cover on ground in forested areas. 

All selected land cover classes from the LC2000 dataset are initially extracted as binary masks for each class at the native pixel resolution. Certain land cover classes are subsequently weighted by area according to the percentage of tree cover, as provided by the C3S Land Cover classification. The weighted LC2000 classes are then spatially aggregated to match the grid resolution of the respective Snow Cover Extent (SCE) products: 0.01° × 0.01° for MODIS and SLSTR, and 0.05° × 0.05° for AVHRR-based SCE retrievals. The aggregated classes are combined to form a Land Cover Density (LCD) layer.

Table 2.2: Classes and original codes of the LC2000 data set used for the preparation of the transmissivity map and applied weighting factor.

 Surface Class Maps

For the pre-classification of snow free areas using a Normalized Difference Snow Index Threshold map (see Section 2.2.3), three surface class maps (SCM) were generated based on the LC2000 dataset. The following LC2000 classes were prepared in the same way as the water and permanent snow and ice area mask (see Section 2.2.1.1):

Table 2.3: Classes and original codes of the LC2000 data set used for the preparation of the NDSI threshold map (see Section 2.2.3).

The SCM are prepared as follows:

\begin{align}
&SCM1 = A + B \tag{Eq. 2.1} \\
&SCM2 = C + D + E + F  \tag{Eq. 2.2} \\
&SCM3 = G  \tag{Eq. 2.3} 
\end{align}

Digital Elevation Model

A global Digital Elevation Model (DEM) with a spatial resolution of 90 meters, derived from TanDEM-X data acquired between 2011 and 2015, was released by the German Aerospace Center (DLR) in 2018. However, detailed evaluation of the publicly available version revealed that it is neither void-filled nor subject to comprehensive quality control. Consequently, for regions between 60°N and 60°S, the 90-meter resolution DEM from the NASA Shuttle Radar Topography Mission (SRTM), version 4.1 (available at http://srtm.csi.cgiar.org), based on data acquired in 2000, was selected as the primary elevation source (Reuter et al., 2007).

For land areas located north of 60°N and south of 60°S, ellipsoidal heights from the TanDEM-X DEM were converted to orthometric heights to ensure consistency with SRTM v4.1 and ASTER GDEM v2, both of which reference the WGS84/EGM96 geoid. Following this conversion, void-filling procedures were applied to the TanDEM-X orthometric height data.

Due to the presence of spurious elevation gradients in the TanDEM-X DEM over large inland water bodies in northern latitudes, additional corrections were implemented. Water bodies larger than approximately 10 km² north of 60°N were identified using the C3S Land Cover classification. In these identified areas, TanDEM-X elevations were replaced with data from the ASTER GDEM v2, a global DEM product jointly developed by NASA and METI, to improve elevation accuracy over water surfaces.

Normalized Difference Snow Index Threshold Map

A Normalized Difference Snow Index (NDSI) threshold map is used for the pre-classification of snow free areas. The concept is based on the premise that snow occurrence is unlikely in regions closer to the equator, at lower elevations and on certain land cover types.

The NDSI threshold map incorporates

1) the geographic latitude Φ ,  

NDSIthr_1(\Lambda,\phi) = \left\{%   
\begin{array}{ll}
    NDSI_{min}, & \phi_1 \leq \phi \leq \phi_0; \\
    m_1 * \phi + b_1, &  \phi_2 < {\phi} < \phi_1; \\
    NDSI_{max}, & \phi_3 \leq \phi \leq \phi_2; \\
    m_2 * \phi + b_2, & \phi_4 < {\phi} < \phi_3; \\
	NDSI_{min}, & \phi_5 \leq \phi \leq \phi_4; \\
\end{array} \tag{Eq. 2.4}  %   
\right.

with

\begin{align}
&\phi = [-90, 90]; \Lambda = [-180, 180]; NDSI_{min} = -0.10;  NDSI_{max} = +0.40; b_1 = 54; b_2 = -54; \\
&\phi_0 = 90;  \phi_1 = 58; \phi_2 = 38; \phi_3 = -38; \phi_4 = -58; \phi_5 = -90;
\end{align}

and

\begin{align}
m_1 &= \frac{\phi_2 - \phi_1}{NDSI_{max} - NDSI_{min}}  \tag{Eq. 2.5} \\
m_2 &= \frac{\phi_4 - \phi_3}{NDSI_{min} - NDSI_{max}}  \tag{Eq. 2.6}
\end{align}

2) the terrain elevation z derived from the DEM ( Section 2.2.2 ),

NDSIthr_2(z) = \left\{%
\begin{array}{ll}
    &NDSIthr_1(\Lambda, \phi), &  z \leq z_0; \\
	&NDSIthr_1(\Lambda, \phi) - (z - z_0) * 20E^{-5}, & z > z_0; \\
\end{array}  \tag{Eq. 2.7}  %
\right.

with

z_0 = 500;

 resulting values smaller than

NDSI_{min}

are reset:

NDSIthr_2 = \left\{%
\begin{array}{ll}
    NDSIthr_2(z), &  NDSIthr_2=[-0.10, 0.40]; \\
	-0.10, & NDSIthr_2 < -0.10; \\
\end{array}  \tag{Eq. 2.8}  %
\right.

3) surface class maps based on the LC2000 dataset (see Section 2.2.1.4 ).

NDSIthr_3(SCM1) = NDSIthr_2 + \frac{SCM1}{100} * 0.20  \tag{Eq. 2.9}  

with

SCM1 = [0, 100]

Y(\phi) = \left\{%
\begin{array}{ll}
    y_{min}, & \phi_1 \leq \phi \leq \phi_0; \\
	m_3 * \phi + b_3, & \phi_2 < {\phi} < \phi_1; \\
    y_{max}, & \phi_3 \leq \phi \leq \phi_2; \\
	m_4 * \phi + b_4, & \phi_4 < {\phi} < \phi_3; \\
	y_{min}, & \phi_5 \leq \phi \leq \phi_4; \\
\end{array}  \tag{Eq. 2.10}  %
\right.

with

\begin{align}
&y_{min} = 0.00;
&y_{max} = 0.20;\\
&b_3 = \phi_1;
&b_4 = \phi_4;\\
&m_3 = \frac{|\phi_2 - \phi_1|}{y_{max} - y_{min}} \tag{Eq. 2.11} \\ 
&m_4 = \frac{|\phi_4 - \phi_3|}{y_{min} - y_{max}} \tag{Eq. 2.12}
\end{align}

NDSIthr_3(SCM2) = NDSIthr_3(SCM1) + \frac{SCM2}{100} * Y(\phi);  \tag{Eq. 2.13}   

with

SCM2 = [0, 100].

NDSIthr = \left\{%
\begin{array}{ll}
	-0.10, & SCM3 > 0; \\
	NDSIthr_3(SCM2), & SCM3 = 0;
\end{array}  \tag{Eq. 2.14}   %
\right.

with

NDSIthr = [-0.10, 0.60].

The basic NDSI threshold map representing winter conditions on both hemispheres is shown in Figure 2.1.

Figure 2.1: Basic NDSI threshold map.

Two-way Forest Transmissivity

To retrieve in forested areas the snow cover fraction on ground, a canopy correction layer representing the transmissivity of the forest canopy is generated. This layer is derived from (i) tree cover density (TCD) data obtained from Landsat imagery acquired in the year 2000, with a spatial resolution of 30 meters, published by Hansen et al., 2013, available at https://data.globalforestwatch.org/datasets/14228e6347c44f5691572169e9e107ad and (ii) the LCD based on aggregated forest classes from the LC2000 data set (see Section 2.2.1).

The estimation of the forest two-way transmissivity is based on an asymmetric sigmoidal fit function, based on the combined TCD and the LCD layers:

t^2 (f)= a + \frac{ b - a}{(1+\frac{f}{c}^d )^e}  \tag{Eq. 2.15}   

with

f = TCD * LCD  \tag{Eq. 2.16}   

and the parameters a = -0.250493, b = 0.9836593, c = 158975900.0, d = 0.5359928, and e = 2898.161.

The resulting domain range 

D[t_{min}^2, t_{max}^2]

is linearly stretched to the co-domain W = [0.06, 1] for AVHRR, and to the co-domain W = [0.08, 1] for MODIS and SLSTR:

t^2(x,y) = 1-\frac{((1-W_{min})*(t_{max}^2- t^2))}{t_{max}^2- t_{min}^2}  \tag{Eq. 2.17}   

The variance of the transmissivity per pixel, used for the uncertainty estimation per observed land pixel (see Section 3.4), is estimated by

S^2_{t^2}= a * (t^2)^2 + b * t^2 + c  \tag{Eq. 2.18}   

with

a = 6.1E^{-3}, b = 5E^{-4},

and

 c = 6E^{-5}.

Spectral Reflectance Maps

Spectral reflectance maps are used as input for the SCE retrieval (Section 3). The maps represent the spectral characteristics of snow free ground

(ρ_{\lambda,ground})

and snow free forest

(ρ_{\lambda,forest})

for a given sensor during the main winter season. The preparation is based on statistical analysis of the spectral reflectance of the first 30 snow free satellite observations per pixel after the winter season per year. The period for extracting snow free spectral reflectance values per pixel starts on January, 1st per year for the Northern Hemisphere, and on July, 1st per year for the Southern Hemisphere for the years listed in Table 2.4.

Table 2.4: Selected years per satellite sensor used for the generation of spectral reflectance maps.

For n snow free cloud free observations per year per pixel, the quartiles (k = 1, 2, 3) of the spectral reflectance are defined for every land pixel:

Q_k = ((k * (n + 1)) / 4), \hbox{with } n = [0, 30]  \tag{Eq. 2.19}   

Outliers are excluded from the spectral reflectances used for the further procedure:

\rho_{\lambda,obs}(x,y,t) = (Q_1 - IQR * 1.5, Q_3 + IQR * 1.5)  \tag{Eq. 2.20}   

with

IQR = Q_3 - Q_1.  \tag{Eq. 2.21}   

Non-forested areas

For non-forested areas, defined by

t^2 = 1

, the spectral reflectance map is directly taken from the observed spectral reflectances per sensor at snow free conditions for every selected year:

\rho_{\lambda,obs\_ground}(x,y) = min(\frac{1}{n}\sum\limits_{i=Q_1}^{Q_2}\rho_{\lambda,obs}(x,y,t_i)).  \tag{Eq. 2.22}   

For some pixels, the pre-condition of 30 snow free and cloud free observations cannot be fulfilled. To cover these pixels as well in the spectral reflectance map representing snow free ground, the spectral reflectance of the observed pixels 

\rho_{\lambda,obs\_ground}(x,y)

are gradually extrapolated, using a hierachical filter window size of 

[k \hbox{ x } k]

 pixels, with

k \in \{5, 11, 21, 31, 41, 51, 61, 71, 81, 91, 141\}.

The extrapolated values

\rho_{\lambda,ground}(x,y)

are used for any missing non-forested land pixels in 

\rho_{\lambda,obs\_ground}(x,y).

 Forested areas

In forested areas, the spectral reflectance maps are retrieved based on the assumption that the observed signal is a mixture of reflectance from both the forest canopy and the underlying ground, with their relative contributions determined by the canopy density, represented by the two-way forest transmissivity (see Section 2.2.4).

In forested areas, the initial spectral reflectance is obtained from the observed spectral reflectances, following the exclusion of outlier values: 

\rho_{\lambda,obs\_forest}(x,y) = min(\frac{1}{n}\sum\limits_{i=min}^{Q_1}\rho_{\lambda,obs}(x,y,t_i)),  \tag{Eq. 2.23}   

For pixels with a tree cover density greater than 50%, represented by 

t^2 < 0.23

the observed spectral reflectance 

\rho_{\lambda,obs\_forest}(x,y)

is extrapolated and tailored to all forested pixels 

\rho_{\lambda,forest}(x,y) = \rho_{\lambda,obs\_forest}(x,y, t^2 < 1)  \tag{Eq. 2.24}   

using the same filter window sizes as described in Section 2.2.5.1.

For the signal contribution from the ground underneath the forest canopy, the extrapolated spectral reflectance map for snow free ground

\rho_{\lambda,ground}(x,y)

is used as initial spectral reflectance.

The extrapolated spectral reflectance components for the forest canopy

\rho_{\lambda,forest}(x,y)

and underlying ground

\rho_{\lambda,ground}(x,y)

are iteratively combined as a function of the two-way transmissivity

t^2(x,y)

until the residual between the assimilated and the observed spectral reflectance falls within an acceptable margin of error:

\Delta \rho_{\lambda}(x,y) = \bigg(min\big(\frac{1}{n}\sum\limits_{i=Q_1}^{Q_2}\rho_{\lambda,obs}(x,y,t_i)\big)\bigg) - \bigg(t^2(x,y) * \rho_{\lambda,ground, iter}(x,y) +  \big(1 - t^2(x,y)\big) * \rho_{\lambda,forest, iter}(x,y)\bigg) \approx 0  \tag{Eq. 2.25}   

with

iter = [0, 150]

or

\Delta \rho_{\lambda}(x,y) < 1E^{-9}.

In forested areas, the iteratively derived spectral reflectance maps for snow-free ground and snow-free forest along with the two-way forest transmissivity map are utilized as auxiliary data set for the fractional snow cover classification per pixel (see Section 3.4).

The associated variances of these spectral reflectance maps are estimated per sensor (Table 2.5).

Table 2.5: Variances of spectral reflectance maps per sensor.

Algorithms 

This section provides a detailed description of the algorithms used for the SCE retrieval, including the associated uncertainty characterization. An overview on the processing chain is provided in Figure 3.1. The approach contains five majore modules, which are described in the following sub-sections.

Figure 3.1: Processing chain for the Snow Cover Fraction product generation, including viewable snow and snow on ground corrected for the shading effect of the forest canopy.

Module 1 - Preprocessing of satellite data 

 AVHRR

Daily global NOAA AVHRR spectral at-satellite radiance data and solar and sensor zenith and azimuth angles in degrees are imported with 0.05 x 0.05 deg. pixel spacing. The reflectance data are converted using the associated scale and offset factors and corrected by division by the cosine of the solar zenith angle. Resulting values are scaled between 0 and 1. For emissive bands, the data are converted into brightness temperatures in Kelvin by applying the associated scale and offset factors.

For channel 3b, the reflective part of the spectral range is extracted as described e.g. by Baum and Trepte (1999) using the parameterisation of Trishchenko (2006).

 MODIS

The following Terra MODIS Collection 6.1 data sets are used:

• MOD021KM: Level 1B Calibrated Radiances - 1 km, swath data provided in 5 minute granules
• MOD03: Level-1 Geolocation product, swath data provided in 5 minute granules

From the 36 discrete spectral bands available on the Terra MODIS instrument, the reflective bands 4 (0.555 µm) and 6 (1.6 µm) and the emissive bands 20 (3.7 µm), 31 (11 µm) and 32 (12 µm) from the MOD021KM data set are extracted and geolocated onto a geographic latitude-longitude grid based on the WGS84 reference ellipsoid, using the geolocation information provided in the MOD03 data set.

Radiances from the reflective bands are converted into unitless top of atmosphere reflectance, while radiances from the emissive bands are converted into brightness temperatures in Kelvin. All data are projected to the same grid with a pixel spacing of 0.01 deg x 0.01 deg. In addition, solar and sensor zenith and azimuth angles in degrees are extracted from the MOD03 data set and geolocated into the same map projection and grid spacing. 

  SLSTR

Extracted Sentinel-3 A/B SLSTR L1B data are provided in frames, including calibrated at-satellite reflectances with about 500 m grid spacing, brightness temperatures with about 1 km grid spacing, and associated auxiliary data with about 16 km grid spacing. 

The extracted Sentinel-3A/B SLSTR Level-1B data are imported and geolocated onto a geographic latitude–longitude grid based on the WGS84 reference ellipsoid, utilizing the corresponding per-frame metadata.

The data of the reflective bands S1 (0.555 µm) and S5 (1.6 µm) are imported and converted into unitless top of atmosphere reflectance data with 0.01 deg. x 0.01 deg. pixel spacing, scaled between the values 0 and 1. The data of the emissive bands S7 (3.7 µm), S8 (11 µm) and S9 (12 µm) are imported and converted into brightness temperatures in Kelvin with a pixel spacing of 0.01 deg. x 0.01 deg. The solar and sensor zenith and azimuth angles in degrees are read and geolocated into the same map projection with a pixel spacing of 0.16 deg. x 0.16 deg.

Module 2 - Cloud screening

The cloud screening is based on different approaches, based on the spectral bands and available cloud products. Clouds over water are masked as water in the final product. All cloud free land pixels are considered for the next processing step, the pre-classification of snow free areas. Clouds are masked in the Snow Cover Fraction products and in the associated uncertainty layer.

AVHRR

The CLARA A3 cloud probability mask is used to mask clouds from NOAA AVHRR data. Clouds are masked if the cloud probability of a pixel is greater than 50%.

 MODIS and SLSTR

The cloud screening from Terra MODIS and Sentinel-3 SLSTR data is based on an adapted version of the Simple Cloud Detection Algorithm (SCDA) (Metsämäki et al., 2015). The decision tree to mask clouds is shown in Figure 3.2.

Figure 3.2: Processing chain of the cloud screening approach, based on the Simple Cloud Detection Algorithm (Metsämäki et al., 2015).

Module 3 - Binary pre-classification of snow free areas 

The pre-classification of snow free pixel is based on thresholds applied on the NDSI per scene and on a brightness temperature band (11 µm) (BT11). The NDSI is calculated based on spectral top of atmosphere reflectance ρ in the visible (vis) and in the short wave (swir) infrared spectral range:

NDSI = \frac{\rho_{vis} - \rho_{swir}}{\rho_{vis} + \rho_{swir}} \tag{Eq. 3.1}

A pixel is classified as snow free, if at least one of the following two rules is fulfilled:

NDSI < NDSIthr \tag{Eq. 3.2}

or

BT11 > BTthr \tag{Eq. 3.3}

with NDSIthr being the threshold for the NDSI value, and BTthr being the threshold for the brightness temperature in Kelvin for a spectral band centered around 11 µm.

The NDSIthr per pixel is based on the NDSI threshold map (Section 2.2.3), temporally adapted by the factors provided in Table 3.1. 

Table 3.1: Temporal variation of NDSI threshold map (NDSIthrmap) (see Section 2.2.3). NH = Northern Hemisphere. SH = Southern Hemisphere.

Sensor specific information of spectral bands used for the calculation of the NDSI and the setting of the Brightness Temperature threshold (BTthr) per sensor are provided in Table 3.2. 

 Table 3.2: Sensor-specific pectral bands and thresholds used for pre-classification of snow free areas. 

\rho_{vis}

\rho_{swir}

bt_{tir}

* reflective part is extracted (see Section 3.1.1).

If a pixel is pre-classified as snow-free, the associated uncertainty is set to zero in the final product.

All remaining pixels are considered for the snow cover fraction algorithm (see Section 3.4).

Module 4 - Fractional Snow Cover Classification

The algorithm of Metsämäki et al. (2012, 2015) is selected as baseline algorithm for retrieving the Snow Cover Fraction (SCF) at the pixel-level from optical satellite data. The SCF classification relies on the generally higher reflectance of snow in the visible spectral range relative to other land cover types. The retrieval is based on the assumption that the observed spectral reflectance is a function of the SCF:

\rho_{obs}(SCF) = (1 - t^2) * \rho_{forest} + t^2 * \big( SCF * \rho_{snow} + (1 - SCF ) * \rho_{ground} \big) \tag{Eq. 3.4}

with

\begin{align}
\rho_{\lambda,obs}&=\hbox{observed spectral reflectance from unit area}\\
SCF&=\hbox{Snow Cover Fraction per unit area in percentage}\\
t^2&=\hbox{apparent forest canopy two-way transmissivity for the unit area}\\
\rho_{\lambda,ground}&=\hbox{spectral reflectance of snow free ground}\\
\rho_{\lambda,forest}&=\hbox{spectral reflectance of snow free forest}\\
\rho_{\lambda,snow}&=\hbox{spectral reflectance of melting snow}\\
\end{align}

This approach exploits the top of atmosphere reflectance (ρ, TOAR) of a visible band, preferably centred around 0.550 µm. The auxiliary spectral reflectance maps ρ _λ,ground   and  ρ _λ,forest  are generated from multi-temporal satellite observations under snow-free conditions, using per-pixel time series data ( Section 2.2.5 ). In regions with recurrent seasonal snow cover, such cloud-free, snow-free acquisitions are primarily available during the melt and summer seasons, when solar illumination is optimal for reliable spectral reflectance retrieval.

To estimate the spectral reflectance of melting snow, 

\rho_{\lambda,snow},

13 Landsat scenes with varying snow conditions in non-forested areas are selected. From these scenes, the SCF is retrieved using a multi-spectral unmixing algorithm with locally adaptive end-member selection (Keuris et al. 2023) and aggregated to match the pixel spacing of the global SCF product. Landsat-based SCF values are then used to solve the inverted SCFG retrieval equation for

\rho_{\lambda,snow}.

The analysis focuses on pixels with high snow cover values (85% ≤ SCF ≤ 100%). Wet or refrozen snow, characterised by large grain sizes, typically exhibits significantly lower spectral reflectances than fresh snow (e.g. Nakamura 2001). From the resulting distribution of snow pixel reflectances, the tenth percentile is defined to represent the spectral reflectance of typical wet snow conditions under full snow cover, yielding: 

\rho_{\lambda,snow} = 0.55. \tag{Eq. 3.5}

A globally estimated variance 

S^2_{snow}(\rho_{\lambda}) = 0.056 \tag{Eq. 3.6}

is used for uncertainty estimation.

To enable the use of the spectrally dependent auxiliary data also during the primary winter season, a correction is applied to the observed reflectance values based on the per-pixel solar zenith angle (θ). This adjustment compensates for the reduced solar illumination at high solar zenith angles:

\rho_{\lambda,obs} =  \left\{%
\begin{array}{ll}\\
\rho_{\lambda,obs} * \cos (\theta - \theta_0) & \theta \geq \theta_0 \\
\rho_{\lambda,obs} & \theta < \theta_0\\
\end{array}%
\right. \tag{Eq. 3.7}

with 

\theta_0 = 50°.

The original SCF retrieval method introduced by Metsämäki et al. (2012, 2015) utilizes a spectrally dependent canopy transmissivity map in combination with globally constant spectral reflectance values for snow-free ground, snow-free dense forest, and wet snow. To address the considerable variability in spectral reflectance across different geographic regions and sensor types, this methodology has been revised. In the adapted approach, spatially explicit spectral reflectance maps are employed to represent snow-free ground and forest canopy conditions at the pixel level for each sensor (see Section 2.2.5). Furthermore, forest canopy correction is performed using a transmissivity map that is independent of spectral reflectance properties. This modified framework allows for consistent retrieval of both the Snow Cover Fraction on Ground (SCFG) and the Snow Cover Fraction Viewable from Above (SCFV). In areas without forest cover, SCFG and SCFV are equivalent.

The Snow Cover Fraction on Ground (SCFG) per pixel accounting for the shading effect of the forest canopy adapted from Metsämäki et al. (2012, 2015) is estimated by:

SCFG =  \frac{\frac{1}{t^2} * ρ_{λ,obs} + (1 - \frac{1}{t^2})*ρ_{λ,forest}- ρ_{λ,ground}}{ρ_{λ,snow} - ρ_{λ,ground}} \tag{Eq. 3.8}

To calculate the Snow Cover Fraction Viewable (SCFV) from above by the satellite sensor, the transmissivity is set to 

t^2 = 1

for all land pixels. Thus, the equation for SCFG retrieval can be simplified for SCFV to:

SCFV =  \frac{ ρ_{λ,obs} - ρ_{λ,background} }{ ρ_{λ,snow}- ρ_{λ,background} } \tag{Eq. 3.9}

with

ρ_{λ,background} = (1 - t^2 ) * ρ_{λ,forest} + t^2 * ρ_{λ,ground} \tag{Eq. 3.10}

The Snow Cover Fraction (SCF) uncertainty for both, SCFG and SCFV, is quantified as the unbiased Root Mean Square Error (RMSE) for each observed land pixel. The RMSE for a general measurand  y can be decomposed into variance and bias components (Scharf, 1991):

RMSE(y) = \sqrt{MSE(y)} = \sqrt{variance(y) + bias(y, y_{ref})^2} \tag{Eq. 3.11}

with y_ref being a reference data set, usually consisting of an extended sample of ground truth data.

For the global SCE CDR, the availability of ground truth data that adequately represents the range of climatic regions and the full temporal coverage is limited. Continuous in-situ snow measurements exist at only a few locations worldwide, with the majority concentrated in the Northern Hemisphere. To avoid introducing additional uncertainties and regional biases, a global bias of zero is assumed. Consequently, SCF uncertainty is expressed as the unbiased RMSE, i.e.

RMSE_{unbiased}(SCF) = \sqrt{variance(SCF)} = \sqrt{S^2_{SCF}} \tag{Eq. 3.12}

The uncertainty characterisation for the SCF estimation is based on the law of propagation of uncertainties (e.g. Taylor, 1982). The variance S² of the SCF estimate can be described as the statistical uncertainty of the SCF (E_stat)² and is the sum of uncertainty contributions from the input parameters per pixel, following the approach described by Metsämäki et al (2015) and Salminen et al. (2018):

S^2_{SCF}(\rho_{\lambda,obs}(SCF), t^2) = (E_{stat})^2 = E^2_{obs} + E^2_{t^2} + E^2_{snow} + E^2_{forest} + E^2_{ground} \tag{Eq. 3.13}

The uncertainty contribution from the actual spectral observation includes any uncertainty caused by the input satellite data, such as uncertainties in the geolocation, radiometric uncertainty, impact of the atmosphere on the radiance measured at the satellite, or the correction of the spectral reflectance in dependence on the solar zenith angle. The uncertainty contribution of the observation is assumed to be covered by the other uncertainties (Salminen et al., 2013, 2018), and is assumed to be zero in the statistical SCF uncertainty characteristation.

The uncertainty contribution from the two-way transmissivity is retrieved by

\begin{align}E^2_{t^2} = \big(\frac{\partial(SCF)}{\partial t^2} \big)^2 S^2_{t^2} 
= \big(\frac{1}{t^2}\frac{(\rho_{forest} - \rho_{ground}) + SCF(\rho_{ground} - \rho_{snow})}{\rho_{snow} - \rho_{ground}}\big)^2 S^2_{t^2}
\end{align} \tag{Eq. 3.14}

using the variance of the transmissivity described in Section 2.2.4.

The uncertainty contribution from the melting snow parameter is obtained by

E^2_{snow} = \big(\frac{\partial(SCF)}{\partial \rho_{snow}} \big)^2 S^2_{snow} = \big(\frac{SCF}{\rho_{snow} - \rho_{ground}}\big)^2 S^2_{snow} \tag{Eq. 3.15}

using the globally estimated variance for melting snow as described above.

The uncertainty contribution from the snow-free forest reflectance is retrieved by

E^2_{forest} = \big(\frac{\partial(SCF)}{\partial \rho_{forest}} \big)^2 S^2_{forest} = \big(\frac{1 - \frac{1}{t^2}}{\rho_{snow} - \rho_{ground}} \big)^2 S^2_{forest} \tag{Eq. 3.16}

using the globally estimated variance for the spectral reflectance of snow-free forest per sensor (Table 2.5).

The uncertainty contribution from the snow-free ground reflectance is retrieved by

E^2_{ground} = \big(\frac{\partial(SCF)}{\partial \rho_{ground}} \big)^2 S^2_{ground} = \big(\frac{SCF - 1}{\rho_{snow} - \rho_{ground}} \big)^2 S^2_{ground} \tag{Eq. 3.17}

using the globally estimated variance for the spectral reflectance of snow-free ground per sensor (Table 2.5).

 Module 5 - Post-processing

AVHRR

For AVHRR, the following post-classification procedure is applied to reset mis-classified snow covered pixels resulting from Module 4 (Section 3.4) to snow free:

1. In the tropics, based geographic latitude Φ in degrees, elevation z in meters from DEM (Section 2.2.2), observed unitless spectral top of atmosphere reflectance and brightness temperature in Kelvin per pixel: 

   SCF = \left\{% \begin{array}{ll}   0, & \hbox{if } (−15 < \phi < 15 \land z < 1000 ) \land \\ & (\rho_{Ch1,obs}=0.3 \lor Ch5 > 270.0 ) \\   SCFG ∨ SCFV, & \hbox{else} \end{array}% \right. \tag{Eq. 3.18}

2. On Southern Hemisphere, for unreliable spectral top of atmosphere reflectance values: 

   SCF = \left\{% \begin{array}{ll}   0, & \hbox{if }\phi < 0 \land \\ & \big((\rho_{Ch1,obs} - \rho_{Ch3b\_refl,obs}) \leq 0.2 \land \rho_{Ch3b\_refl,obs} > 0.1 \big) \\   SCFG ∨ SCFV, & \hbox{else} \end{array}% \right. \tag{Eq. 3.19}

3. Global, for unreliable spectral top of atmosphere reflectance values: 

   SCF = \left\{% \begin{array}{ll}   0, & \hbox{if } \big((\rho_{Ch1,obs} - \rho_{Ch3b\_refl,obs}) < 0.1 \land \rho_{Ch3b\_refl,obs} > 0.2 \big) \\ SCFG ∨ SCFV, & \hbox{else} \end{array}% \right. \tag{Eq. 3.20}

 MODIS and SLSTR

The MODIS and SLSTR based SCF products are generated per granule and frame, respectively. The term "frame" is used in the following a synonym for both, granule and frame. The SCF processing chain is applied on each frame. For each pixel, the solar zenith angle and sensor zenith angle are checked. Pixels with solar zenith angles (sza) greater than 83° are classified as night or polar night. For pixels with a sensor zenith angle (vza) greater than 65° and a solar zenith angle lower than 83, a reliable classification is not possible. Such pixels are thus classified as retrieval / classification failed. In case, no input data are available for a pixel, these are flagged as no satellite acquisition.

When all individual frames for one day have been processed, the SCF tiles are merged to one global map per day. For overlapping pixels in neighboured frames, a set of hierarchical merging criteria (Table 3.3) are applied to set the output pixel classification. The same merging criteria are used for the SCFG and the SCFV products, as well as for the associated uncertainty layers.

Table 3.3: Merging criteria for SCF classifications from overlapping frames of the same date.

*) If a pixel is observed in both frames, the observation with the smaller vza is used if the difference between the observed sza is lower than 20° and the difference between the observed vza is lower than 40°.

 Snow to cloud conversion

The cloud mask in the daily product is updated based on mis-classified snow covered pixels resulting from Module 4 (Section 3.4) considering for each pixel the last observed classification, the 2 meters temperature and preciptiation data from ERA5Land data downscaled to the grid spacing of the SCF products averaged over the satellite acquisition times. Snow can only fall if the meteorological data since the satellite acquisition time of the previous day d_-1 and of the current day d₀ fulfill the physical pre-conditions:

SNOW_{met} = \left\{%
\begin{array}{ll}
	0, & \big(\bar{T}(d_{-1}, d_0) > 273.15 K \lor \bar{P}(d_{-1}, d_0) < 0.003 m \big)\\
	1, & \big(\bar{T}(d_{-1}, d_0) \leq 273.15 K \lor \bar{P}(d_{-1}, d_0) \geq 0.003 m \big)
\end{array}%
\right. \tag{Eq. 3.21}

In the following cases, a pixel classified on the current date (d₀) as snow is reset to cloud in the SCF product (SCF(x,y,d₀) and in the associated uncertainty (SCF_unc(x,y,d₀) of the pixel:

SCF(x,y,d_0) \land SCF_{unc}(x,y,d_0) = \hbox{Cloud}\quad\left\{%
\begin{array}{ll}
	& SCF(x,y,d_{-1}) = 0 \land SCF(x,y,d_0) > 0 \land SNOW_{met} = 0 \\
	& SCF(x,y,d_{-1}) = Cloud \land SCF(x,y,d_0) > 0 \land SNOW_{met} = 0 \\
    & SCF(x,y,d_0) > 0 \land \bar{T}(d_{-1}, d_0) > 298.15 K \\
\end{array}%
\right. \tag{Eq. 3.22}

 Overlay of static masks

Static masks for water bodies, permanent snow and ice areas (Section 2.2.1.1) and salt lakes (Section 2.2.1.2) are overlaid on the daily SCFG and SCFV products and the associated uncertainty layers.

Output data

The daily global SCF products contain the snow cover fraction per pixel, and an associated uncertainty layer providing the unbiased RMSE per observed pixel. SCF products before 24 February 2000 are based on NOAA AVHRR/2 data with a pixel spacing of 0.05 deg x 0.05 deg. SCF products from 24 February 2000 onwards are based on Terra MODIS data (until 31 December 2022) and Sentinel-3 A&B SLSTR data (from 1 January 2023 onwards), and have a pixel spacing of 0.01 deg x 0.01 deg. 

The SCFV product file provides in forested areas the information on the snow on the top of the forest canopy, i.e. snow viewable from above.

The SCFG product file provides in forested areas the information on the snow on ground, retrieved by applying a canopy correction.

Pixels classified as clouds, water bodies, salt lakes, permanent snow and ice areas, (polar) night or erroneous pixels are masked in all product layers.

SCFG and SCFV product examples and the associated uncertainty layers are shown in Figure 4.1 and Figure 4.2, respectively.

Figure 4.1: a) SCFG product example and b) associated unbiased Root Mean Square Error from NOAA-7 AVHRR/2 data of 20 March 1982.

 Figure 4.2: a) SCFV product example and b) associated unbiased Root Mean Square Error from Sentinel-3A&B SLSTR data of 20 March 2023.

Daily SCFG and SCFV products contain additionally information on the scanline time per sensor as fractional hour of the day, and on the sensor zenith angle per pixel in degrees. 

A detailed description of the output products is provided in the  Product User Guide and Specifications.

References

Baum, B. A., and Q. Trepte. 1999. A Grouped Threshold Approach for Scene Identification in AVHRR Imagery. Journal of Atmospheric and Oceanic Technology, 16:793–800. doi:https://doi.org/10.1175/1520-0426(1999)016<0793:AGTAFS>2.0.CO;2.

ESA Land Cover CCI project team; Defourny, P. 2019. ESA Land Cover Climate Change Initiative (Land_Cover_cci): Global Land Cover Maps, Version 2.0.7. Centre for Environmental Data Analysis, 2019. https://catalogue.ceda.ac.uk/uuid/b382ebe6679d44b8b0e68ea4ef4b701c

Hansen, M. C., P. V. Potapov, R. Moore, M. Hancher, S. A. Turubanova, A. Tyukavina, D. Thau, S. V. Stehman, S. J. Goetz, T. R. Loveland, A. Kommareddy, A. Egorov, L. Chini, C. O. Justice, and J. R. G. Townshend. 2013. High-resolution global maps of 21st-century forest cover change. Science, 342:850–853. doi:10.1126/science.1244693.

Karlsson, K.-G., M. Stengel, J. F. Meirink, A. Riihelä, J. Trentmann, T. Akkermans, D. Stein, A. Devasthale, S. Eliasson, E. Johansson, N. Håkansson, I. Solodovnik, N. Benas, N. Clerbaux, N. Selbach, M. Schröder, and R. Hollmann. 2023. CLARA-A3: The third edition of the AVHRR-based CM SAF climate data record on clouds, radiation and surface albedo covering the period 1979 to 2023. Earth System Science Data, 15:4901–4926. doi:10.5194/essd-15-4901-2023.

Keuris, L., Hetzenecker, M., Nagler, T., Mölg, N. and Schwaizer, G. 2023. An Adaptive Method for the Estimation of Snow-Covered Fraction with Error propagation for Applications from Local to Global Scales. Remote Sensing, 15, 1231, doi: 10.3390/rs15051231

Metsämäki, S., Mattila, O.-P., Pulliainen, J., Niemi, K., Luojus, K. and Böttcher. K. 2012. “An Optical Reflectance Model-Based Method for Fractional Snow Cover Mapping Applicable to Continental Scale.” Remote Sensing of Environment, 123 (August). Elsevier Inc.: 508–21. doi:10.1016/j.rse.2012.04.010.

Metsämäki, S., Pulliainen, J., Salminen, M., Luojus, K., Wiesmann, A., Solberg R. and Ripper, E. 2015. Introduction to GlobSnow Snow Extent products with considerations for accuracy assessment. Remote Sensing of Environment, 156, 96–108. https://doi.org/10.1016/j.rse.2014.09.018

Nakamura, T., O. Abe, T. Hasegawa, R. Tamura, and T. Ohta. 2001. Spectral reflectance of snow with a known grain-size distribution in successive metamorphism. Cold Regions Science and Technology, 32:13–26. doi:10.1016/S0165-232X(01)00019-2.

Reuter H.I, Nelson, A., Jarvis, A. 2007. An evaluation of void filling interpolation methods for SRTM data, International Journal of Geographic Information Science, 21:9, 983-1008. doi: 10.1080/13658810601169899

Salminen, M., J. Pulliainen, S. Metsämäki, K. Böttcher, and K. Heinilä. 2013. MODIS-derived snow-free ground reflectance statistics of selected Eurasian non-forested land cover types for the application of estimating fractional snow cover. Remote Sensing of Environment, 138:51–64. doi:10.1016/j.rse.2013.07.007.

Salminen, M., J. Pulliainen, S. Metsämäki, J. Ikonen, K. Heinilä, and K. Luojus. 2018. Determination of uncertainty characteristics for the satellite data-based estimation of fractional snow cover. Remote Sensing of Environment, 212:103–113. doi:10.1016/j.rse.2018.04.038.

Scharf, L. L. 1991. Statistical Signal Processing. Detection, Estimation, and Time Series Analysis. Addison-Wesley, Reading Mass.

Schwaizer, G., S. Metsämäki, M. Moisander, K. Luojus, M. Heinrich, S. Wunderle, T. Nagler, J. Lemmetyinen, J. Pulliainen, M. Takala, R. Solberg, P. Venäläinen, L. Zschenderlein, and X. Xiao. 2025. ESA CCI+ Snow ECV: Algorithm Theoretical Basis Document, version 5.0, under review.

Taylor, J. R. 1982. An Introduction to Error Analysis. Mill Valley, CA: University Science Books.

Trishchenko, A. P. 2006. Solar Irradiance and Effective Brightness Temperature for SWIR Channels of AVHRR/NOAA and GOES Imagers. Journal of Atmospheric and Oceanic Technology, 23, 198-210. doi:10.1175/JTECH1850.1.

Related articles