# Assimilation of surface reflectance in snow simulations: Impact on bulk snow variables

J. Revuelto, B. Cluzet, N. Duran, M. Fructus, M. Lafaysse, E. Cosme, M.

Dumont

► To cite this version:

J. Revuelto, B. Cluzet, N. Duran, M. Fructus, M. Lafaysse, et al. Assimilation of surface reflectance in snow simulations: Impact on bulk snow variables. Journal of Hydrology, 2021, 603, pp.126966. 10.1016/j.jhydrol.2021.126966. insu-03466311

HAL Id: insu-03466311

https://insu.hal.science/insu-03466311

Submitted on 6 Dec 2021

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>

# Assimilation of surface reflectance in snow simulations: Impact on bulk snow variables

J. Revuelto<sup>a,b,*</sup>, B. Cluzet<sup>b</sup>, N. Duran<sup>c</sup>, M. Fructus<sup>b</sup>, M. Lafaysse<sup>b</sup>, E. Cosme<sup>d</sup>, M. Dumont<sup>b</sup>

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

## ARTICLE INFO

This manuscript was handled by Marco Borga, Editor-in-Chief, with the assistance of Massimiliano Zappa, Associate Editor

**Keywords:**

- Snowpack modelling
- Snow surface reflectance
- Data assimilation
- Particle filter
- Mountain areas

## ABSTRACT

Data assimilation of snow observations significantly improves the accuracy of snow cover simulations. However, remotely-sensed snowpack observations made in areas of complex topography are typically subject to large error and biases, creating a challenge for data assimilation. To improve the reliability of ensemble snowpack simulations, this study investigated the appropriate conditions for assimilating MODIS-like synthetic surface reflectances. We used a simulation system that included the Particle Filter data assimilation technique. More than 270 ensemble simulations involving assimilation of synthetic observations were conducted in a twin experiment procedure for three snow seasons. These tests were aimed at establishing the spectral combination of MODIS-like reflectances that convey the more information in the assimilation system, rendering the most reliable snowpack simulation, and determining the maximum observation errors that the assimilation system could tolerate. The assimilation of the first seven MODIS-like bands, covering visible and near-infrared wavelengths, provided the best scores compared with any other band combination, and thus are highly recommended for use when possible. The simulation system tolerated a maximum deviation from ground truth of 5% without loss of performance. However, the assimilation of the first seven bands of true MODIS surface of reflectance fails on improving simulation results in rouged mountain areas.

## 1. Introduction

Some of the major natural hazards in mountain areas are directly linked to the snowpack distribution and its evolution over time, and include snow avalanches (Schweizer et al., 2003, 2008) and floods in downstream areas (Gaál et al., 2015). Forecasting snowpack evolution in mountain areas is increasingly based on numerical snowpack models (Lehning et al., 2006; Vionnet et al., 2012; Morin et al., 2020). The marked natural spatial and temporal variability in the snowpack in mountain areas (Scipión et al., 2013; Seidel et al., 2016) necessitates accurate simulation of snowpack processes and accurate meteorological inputs. Hence, snowpack modelling is affected by forcing and model errors that result in discrepancies between the real state and the simulated snowpack (Morin et al., 2020). The accumulation of such discrepancies over time decreases forecasting capabilities, leading to large uncertainties in risk management. Consequently, for robust operational use the snowpack models need to provide upgraded real time forecasts

(Wayand et al., 2015).

Ensemble approaches are commonly used to quantify the uncertainties associated with simulations (e.g., Swinbank et al., 2016; Vernay et al., 2015). Moreover, data assimilation techniques enable optimal combinations of simulations (usually an ensemble representing possible simulation states) and observations to mitigate uncertainties arising from both sources of information (Andreadis and Lettenmaier, 2006; Clark et al., 2008; Largeron et al., 2020). Data assimilation techniques in snowpack modelling are known to significantly reduce simulation uncertainties and biases (Magnusson et al., 2017; Piazzi et al., 2019). However, the possible assimilation observations and their optimal frequencies and spatial extent may vary depending on the assimilation technique, and the characteristics of the observations. In situ measurements are the most accurate observations (Ménard et al., 2019; Winstral et al., 2019), but their generally poor spatial extent and representativeness makes them inappropriate for large scale simulations at reasonable spatial resolution (i.e., several hundreds or thousands of

<footnote><sup>a</sup> Instituto Pirenaico de Ecología, Consejo Superior de Investigaciones Científicas (IPE-CSIC), Zaragoza, Spain</footnote>

<footnote><sup>b</sup> Univ. Grenoble Alpes, Université de Toulouse, Météo-France-CNRS, CNRM, Centre d'Etudes de la Neige, Grenoble, France</footnote>

<footnote><sup>c</sup> Magellium, Unité Observation de la Terre, Ramonville Saint-Agne, France</footnote>

<footnote><sup>d</sup> UGA, CNRS, Institut des Géosciences de l'Environnement (IGE) UMR 5001, Grenoble F-38041, France</footnote>

<footnote>* Corresponding author.
E-mail address: jrevuelto@ipe.csic.es (J. Revuelto).</footnote>
meters). However, different studies have exploited these observations with different data assimilation schemes (Magnusson et al., 2017; Piazzı et al., 2018) Remote sensing observations, including those involving optical space-borne sensors, collect information on the snowpack sur- face at high spatial resolution over large areas (Immerzeel et al., 2009; Gascoin et al., 2015; Dumont and Gascoin, 2016). Although the fre- quency of optical remote sensing in some regions is limited by cloud cover (Hall and Riggs, 2007), the assimilation of sparse satellite obser- vations into snowpack models has significantly improved the simulation results in theoretical experiments based on synthetic observations (Charrois et al., 2016; Cluzet et al., 2020a). Similarly active satellite sensors have been assimilated in large study areas (Cortés et al., 2016; Larue et al., 2018; Margulis et al., 2019), and in some cases combined with in-situ observations (Piazzı et al., 2019), showing encouraging results. These works demonstrate the effort of the snow modeling community on assimilating snow observations to improve forecasting capabilities of snow models. Probably remote sensing techniques capa- bility of retrieving information over extended (and remote) areas is the main reason that has motivated this attempt.

Snow spectral reflectance varies substantially with wavelength in the solar spectrum (Warren, 1982). The information held in the various bands of optical satellite sensors can provide information about various snowpack surface properties (including snow aging, snow microstruc- ture, or light absorbing particle content; Painter et al., 2012; Skiles et al., 2018). In addition, snow cover monitoring using satellites, particularly in highly heterogeneous areas having patchy snow distribution and/or highly complex topography, is affected by a number of sources of un- certainty (Frei et al., 2012; Cluzet et al., 2020a). When this information is assimilated using a data assimilation scheme, the potential deviations of the observation must be taken into account and appropriately rep- resented in the simulation system (Janjic et al., 2018). Various data assimilation techniques have been tested, depending on the character- istics of the snowpack model and the type of observations. The Ensemble Kalman Filter scheme (Evensen, 2003) has been shown to efficiently assimilate remote sensed and ground based snow observations, improving simulation accuracy (Piazzı et al., 2019). However, the implementation of this technique is challenging, and is associated with significant limitations for detailed snowpack models (Piazzı et al., 2018). The Particle Filter (PF) (van Leeuwen, 2009) is a data assimila- tion technique that does not rely on physical assumptions of the system and is well suited to snowpack models (Charrois et al., 2016; Piazzı et al., 2018). The efficiency of the PF in assimilating observations into snowpack models has been tested in many study areas and with differing observation datasets (Thirel et al., 2013; Magnusson et al., 2017). The PF technique largely reduces the uncertainty in ensemble simulations of the snowpack when reproducing the snow depth (SD) and snow water equivalent (SWE) (Smyth et al., 2019), and has also showed encouraging results in distributed simulations (Baba et al., 2018). Moreover, the PF technique only selects (or discards) particles, resampling those members selected depending on the distance to the observed variable, and in this way maintains the consistency between state variables (Magnusson et al., 2017). These attributes have encouraged snow modelling scien- tists to implement the PF technique in assimilating snowpack variables retrieved using satellite sensors (Largeron et al., 2020). Nevertheless, in view of the fast development of both, remote sensing and models arising in the snow science community (Girotto et al., 2020), an assessment on measurement errors that may be admissible for improving simulation results is required.

This study relied on the detailed Crocus snowpack model (Vionnet et al., 2012), which has been used for more than 25 years in the French mountain ranges to provide operational support for avalanche hazard forecasting. To reduce discrepancies between the observed and simu- lated snowpack (Revuelto et al., 2018), recent developments in this simulation system have aimed to assimilate different snowpack vari- ables with PF assimilation algorithm (Charrois et al., 2016; Charrois 2017; Cluzet et al., 2020a,b). The assimilation of synthetic observations,

including SD and snow surface reflectance has already demonstrated a remarkable improvement on the forecasting capabilities of Crocus (Charrois et al., 2016). Similarly, the assimilation of ground based sur- face reflectance also reduced Crocus simulations uncertainties, but oppositely the assimilation of true satellite-retrieved optical reflectance with MODIS sensors failed (Charrois, 2017). This is mainly attributed to the bias of MODIS snow surface reflectances likely due to complex terrain effects (Cluzet et al., 2020a). Previous studies have paved the way for substantial improvements to detailed snowpack simulations, but have also highlighted the need of determining the maximum observa- tion uncertainties that can be tolerated in the assimilation system and which spectral observations convey the more information in the system. This study aims at determining the minimum accuracy of snow surface reflectance satellite observations, needed to guarantee an effective assimilation that improves snow simulation capabilities.

This study describes assimilation experiments using both synthetic and MODIS surface reflectance with CrocOc1.0 (an ensemble data assimilation system; Cluzet et al., 2020b). First, we investigated the performance of PF algorithm with respect to the spectral bands assimi- lated, which enabled us to determine whether there was a particular band combination rendering better simulation scores. Second, we assessed the maximum observation errors acceptable with respect to improving the forecasting capability of simulations using data assimi- lation. To address these objectives we tested the performance of PF in assimilating synthetic observations using a twin experiment involving data assimilation dates selected from true MODIS images made under ideal conditions (no clouds and high viewing zenith angles) for retrieving snow surface information. In view of the potential for implementation of this system to assimilate real optical surface re- flectances, the synthetic observations involved MODIS-like surface re- flectances extracted from various reference members considered to represent the truth snowpack conditions (termed “truth members”). This procedure included information on all snowpack variables at any time, enabling full evaluation of the assimilation, and also overcame potential shortcomings of real observations. The performance of the assimilation system was evaluated by assessing the impact that data assimilation had on the bulk snow variables SD and SWE. Finally, for the band combi- nation rendering the best results, the assimilation of true MODIS surface reflectance was evaluated analyzing the capability of the simulation system on reproducing the observed SD.

The simulation exercise presented here spanned three snow seasons (2014–15, 2015–16, and 2016–17) at the Col du Lautaret study site in the French Alps. In total, 270 data assimilation experiments were carried out, and combined differing MODIS-like spectral bands and various biases and errors introduced in the assimilated reflectances. Section 2 describes the simulation and assimilation systems, provides detailed information about the study period and the study site, and describes the evaluation procedure. Section 3 presents and describes the main results obtained. Section 4 discusses the main results, and Section 5 presents the conclusions from the study.

## 2. Data and methods

The snowpack data assimilation system aims to reduce the dispersion of an ensemble of simulations towards the true evolution of the snow- pack. Globally the simulation system reproduces the snowpack along the snow season. When a snowpack observation is available, the ensemble of simulations is stopped and resampled trough the particle filter. After- wards, the simulation of the resampled ensemble continues until the next assimilation or the end of the snow season if no more observations are available. This section, after describing the study area and period, provides an overview on the simulation system, comprising snowpack simulations, meteorological forcing, generation of the ensemble and the assimilation procedure. Later on, once MODIS data processing is described (Section 2.4), the twin experiment procedure followed in this work is detailed (Section 2.5). Finally the assimilation experiments that
enabled the assessment on the efficiency of the system in view to the different bands combinations assimilated (Section 2.6) and the bias of surface reflectance tolerated by the system (Section 2.7) are specified.

## 2.1. Study site and period

The simulations were conducted in Col du Lautaret, a mountain pass located in the Central French Alps (45.02°N, 6.46°E) between four Alpine massifs (Ecrins, Pelvoux, Grandes Rousses, and Thabor) at an elevation of 2058 m a.s.l. The study area is completely covered by alpine

grassland and has a flat topography and without any predominant exposition (see summer and winter pictures of Fig. 1). Despite the surrounding topography is rugged but locally open (see Fig. 1) the satellite observations of this study site are also influenced by the diffuse radiation scattered by the neighborhood that characterizes snow covered rugged terrain and by direct re-illumination effects (Lamare et al., 2020). Because of its accessibility and high elevation it has been involved in many snow studies (Lamare et al., 2020; Larue et al., 2020; Revuelto et al., 2020; Tuzet et al., 2020). Moreover, this site is equipped with an automatic weather station, acquiring different variables including SD.

<figure><img src="image_4.png" /><figcaption>**Fig. 1.** Col du Lautaret study site location and main topographic characteristics of the surrounding area. This figure, also includes two pictures of the study site in winter (snow covered) and summer (snow free), showing the homogenous topography of the grid cell where simulations were conducted.</figcaption></figure>
As other data assimilation studies in mountain areas (Smyth et al., 2019; Winstral et al., 2019), our study period included three snow seasons (2014–15, 2015–16, and 2016–17), and thus included a wide variety of meteorological events and snowpack conditions.

## 2.2. Open loop simulation system

The simulation platform used for simulating snow evolution is the externalized surface model SURFEX, which integrates the Crocus snowpack model (Vionnet et al., 2012) and the multilayer soil model ISBA-DIF (Decharme et al., 2016). Crocus is a detailed multilayer snowpack model that simulates the most important processes within the snowpack, including all energy and mass exchanges between the snow layers, the soil, and the atmosphere. The TARTES radiative transfer scheme (Libois et al., 2015) was recently implemented in Crocus to provide a vertical profile of solar radiation absorption from snow microstructure profiles (based on the specific surface area: SSA), and the concentrations of light absorbing particles (LAPs; e.g., black carbon and mineral dust) deposited on the snow surface (Tuzet et al., 2017). In this study the multi-physics version of Crocus was used (ESCROC: Ensemble System CROCus; Lafaysse et al., 2017) to account for model errors in ensemble simulations. This way, Crocus simulates the temporal evolution of nearly all variables describing the different snowpack layers (Vionnet et al., 2012) including the seven MODIS like surface reflectance bands, and bulk snow pack variables as the SD and the SWE.

The simulations were conducted on a 250 × 250 m pixel with the topographical attributes (slope, elevation, aspect...) of Col du Lautaret. The meteorological inputs were obtained from the SAFRAN reanalyses (Durand et al., 2009; Vernay et al., 2015). SAFRAN reanalysis provides at hourly time step all atmospheric variables needed to run Crocus model, including air temperature, specific humidity, precipitation phase and rate, direct and diffuse radiation, long wave radiation and wind speed. This analysis scheme has a semi-distributed geometry, providing for discrete units of the Alps of about 1000 km², named “massifs”, and included all meteorological variables needed for running Crocus simulations. We considered the meteorological forcings corresponding to the elevation, aspect, and slope of the study sites for the four SAFRAN massifs (Oisans, Pelvoux, Grandes Rousses, and Thabor), because of the position of our study site at the intersection of these massifs. From each of these four reanalyses we obtained 50 members perturbed using a stochastic approach (Charrois et al., 2016). Using MOCAGE outputs (which is the Chemistry-Transport Model of MétéoFrance (Modèle de Chimie Atmosphérique de Grande Echelle) Josse et al., 2004), an ensemble of LAP deposition fluxes was generated, and stochastically associated with SAFRAN members. The 200 meteorological members obtained were randomly combined with 200 parameterizations of Crocus through the Eltartes option of ESCROC, as described by Cluzet et al. (2020a). An example of the temporal variability of the ensemble of snowpack simulations without data assimilation (open loop), based on the evolution of

SD during the study period, is provided in Fig. 2.

## 2.3. Data assimilation algorithm

The PF technique was implemented using the Sampling Importance Resampling (SIR) algorithm (Gordon et al., 1993; van Leeuwen, 2009). The PF technique estimates the distribution of possible states of a system using a probability density function based on the sampling of the members of an ensemble. Each time that an observation is available, all particles of the ensemble are weighted through the likelihood to the observation, by computing the change in the relative importance of the probability density function. This function also takes into account the inherent error in the observations through a covariance matrix R (Cluzet et al., 2020b). Once the ensemble of possible states (particles) has been weighted, particles are resampled. Those having a negligible weight are discarded, and those having higher weights (closer to the observed value) are replicated. The resampling essentially maintains the total number of particles, but outputs a less dispersive ensemble shifted towards the observed value when properly settled (see Fig. 2 in van Leeuwen, 2009). Further details of this technique and the theory behind it are described by van Leeuwen (2009). The practical implementation of the Crocus model and reflectance assimilation is described in detail by Cluzet et al. (2020b).

## 2.4. MODIS data processing

Many studies have demonstrated the usefulness of MODIS images for snow cover mapping in mountain areas at 500 m spatial resolution (Gascoin et al., 2015; Parajka and Blöschl, 2008). Sub-pixel snow monitoring of the snow cover at 250-m spatial resolution was performed using MODimLab software (Dumont et al., 2012; Sirguey et al., 2009). Multispectral fusion between MOD02HKM and MOD02QKM (Sirguey et al., 2008), enabled this software to generate images at 250 × 250 m spatial resolution to derive various snow-ice products, including the first seven bands of MODIS sensor (respectively 460, 560, 640, 860, 1240, 1640, and 2120 nm central wavelengths) and also cloud cover fraction. These variables were obtained for all available MODIS sensor images from Terra platform along the study period.

MODIS images were used to select data assimilation dates. To ensure the highest quality satellite images, we checked MODIS images for Col du Lautaret to find dates on which the site was fully snow-covered, the sky was clear, and the viewing azimuth angle was high (>60°) (Sirguey et al., 2016). For the 2014–15, 2015–16, and 2016–17 snow seasons 16, 14, and 19 assimilation dates (respectively) met these requirements. For these dates, snow surface reflectance for the first seven MODIS bands was also kept for assimilation.

<figure><img src="image_5.png" /><figcaption>**Fig. 2.** Snow depth evolution of the open loop ensemble for 2014–15 (left), 2015–16 (center), and 2016–17 (right). Two envelopes for each snow season show all members contained in the 25 to 75 percentiles and in the 10 to 90 percentiles. Similarly, the continuous purple line denotes the temporal evolution of the ensemble median. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>
## 2.5. Twin experiment with true and synthetic observations

Twin experiments are typically used to test data assimilation techniques (Matgen et al., 2010; Dubinkina et al., 2011; Browne and van Leeuwen, 2015). In essence, a twin experiment compares two ensembles: an open loop ensemble without data assimilation and a second ensemble assimilating observations. The data assimilation system has been tested with two different observations datasets: MODIS-like surface reflectance (derived from the open loop ensemble) and true MODIS surface reflectance. These two datasets allowed to explore different assimilation scenarios, having a detailed description of the snowpack at any time and also to evaluate the performance of the simulation system when assimilating true satellite observations.

For the distinct assimilation experiments detailed in the next section, synthetic observation values were only obtained for dates in which MODIS observations had high quality (Section 2.4), from one member of the open loop ensemble. This member, named the truth or reference run, was considered to describe the real snowpack state for the given experiment, and contained temporal information on all variables of the simulation system. When assimilating synthetic observations, the reference run was removed from the ensemble.

The assimilation of both true and synthetic observations, had same data assimilation dates which selection was only based on the quality of satellite images, and thus they reproduces true operational assimilation procedures based on the availability of high quality observations providing a realistic temporal repetitiveness of assimilation dates.

## 2.6. Assimilation experiments

One major advantage of the twin experiment framework is that we were able to conduct various tests assimilating distinct reference runs, including true MODIS observations. This enabled assessment of the efficiency of the system as a function of the location of the “true” state within the ensemble. For the three snow seasons we tested the assimilation of seven reference runs (21 in total) randomly selected from the open loop simulation. Table 1 shows the SD and the SWE percentiles (based on mean values during the period when snow was on the ground in at least one member) for the various reference runs in the open loop ensemble. The random selection of several reference runs for each snow season, illustrate how possible snowpack conditions may be simulated and how data assimilation may impact simulation results for different cases. In turn this procedure allows determining which spectral bands convey the more info in the system.

The final objective is the assimilation of true MODIS observations; this way, the covariance matrix of errors (R) of the PF density function had to account for real MODIS sensor errors. In the case study of twin experiments, this is a diagonal matrix having the following errors: 7.1E-4, 4.6E-4, 5.6E-4, 5.6E-4, 2.0E-3, 1.5E-3, and 7.8E-4 (Wright et al., 2014). Errors reported in this later work were determined comparing MODIS observations with spectral albedo field measurements and as far author knows, these are one of the few datasets available on MODIS spectral bands accuracy over snow. Thereby, our

<table><caption>**Table 1**

Percentile of the truth members randomly selected from the open loop ensemble. The classification was based on the SD and SWE bulk variables.</caption><thead><tr><th rowspan="2">Member extraction</th><th colspan="2">2014-15</th><th colspan="2">2015-16</th><th colspan="2">2016-17</th></tr><tr><th>SD</th><th>SWE</th><th>SD</th><th>SWE</th><th>SD</th><th>SWE</th></tr></thead><tbody><tr><td>i</td><td>65</td><td>67</td><td>45</td><td>43</td><td>83</td><td>81</td></tr><tr><td>ii</td><td>77,5</td><td>82</td><td>70</td><td>63,5</td><td>100</td><td>99</td></tr><tr><td>iii</td><td>41,5</td><td>41</td><td>99</td><td>98,5</td><td>47,5</td><td>50</td></tr><tr><td>iv</td><td>3,5</td><td>2,5</td><td>82,5</td><td>64,5</td><td>0,5</td><td>0,5</td></tr><tr><td>v</td><td>95</td><td>97</td><td>61,5</td><td>48</td><td>12,5</td><td>11</td></tr><tr><td>vi</td><td>87</td><td>76</td><td>89,5</td><td>84</td><td>70</td><td>74,5</td></tr><tr><td>vii</td><td>18</td><td>22</td><td>38</td><td>42</td><td>28,5</td><td>26</td></tr></tbody></table>

assimilation experiments used R matrix with elements from deviations reported by Wright et al., 2014. Since some exploratory analysis of Charrois (2017) suggested that an R matrix with elements multiplied by a factor of 5 obtained more reliable ensembles, we also tested the assimilation of synthetic observations with R matrix with such a multiplicative factor (termed here in after F5).

Additionally, the impact of assimilating four combinations of surface reflectance bands was investigated for each of the 21 synthetic observations evaluated. The first of these combinations comprised the assimilation of the seven bands (configuration A). The second combination included the assimilation of the first and fifth MODIS-like reflectances (460 and 1240 nm; configuration B). The third combination tested (configuration C) aimed to evaluate the impact of assimilating visible and the first very near-infrared band together (460, 560, 640, 860 nm); this involved the first four MODIS-like reflectances. The fourth combination (configuration D) assimilated three near-infrared bands (1240, 1640, 2120 nm). All band combinations and configuration names used here after are detailed in Table 2. These band combinations were aimed at evaluating the performance of the system for several spectral band configurations: spectral reflectances more sensitive to impurities on the snow surface (first four MODIS-like bands); reflectances accounting for snow aging (the three higher wavelength MODIS-like bands; Dozier et al., 2009; Skiles et al., 2018); or those assimilating other band combinations having better scores. The band combination showing the best results has been evaluated with true MODIS surface reflectance, to analyze the reliability of true data assimilation of snow surface reflectance. Whereas the assimilation experiments with synthetic observations were evaluated with SD and SWE forecast of the “true” state (reference run); the assimilation of true MODIS surface reflectance was evaluated with the SD observed at the automatic weather station located within the study area.

## 2.7. Biased surface reflectances

To investigate the robustness of data assimilation to random errors and systematic biases in optical satellite observations, we also ran the following experiments. First, three white noises (±2%, ±5%, and ±10%) having random values within these intervals for the various assimilation dates were introduced to the synthetic surface reflectances. Secondly, the synthetic observations were systematically biased for all assimilation dates using +2%, -2%, +5%, -5%, +10%, and -10% of the truth member values. These experiments were only applied to the best band combination identified by the experiments described in Section 2.6.

## 2.8. Evaluation metrics

The evaluation metrics must assess the performance and quality of

<table><caption>**Table 2**

Names of the configurations tested and their corresponding MODIS band combinations and R matrix used in the assimilation algorithm.</caption><thead><tr><th>Configuration name</th><th>MODIS band combination (central wavelength) [nm]</th><th>R matrix</th></tr></thead><tbody><tr><td>WA (reference run)</td><td>Without Assimilation</td><td>Without Assimilation</td></tr><tr><td>A</td><td>460, 560, 640, 860, 1240, 1640, 2120</td><td>Wright et al., 2014</td></tr><tr><td>A-F5</td><td>460, 560, 640, 860, 1240, 1640, 2120</td><td>Wright et al., 2014</td></tr><tr><td></td><td></td><td>x5</td></tr><tr><td>B</td><td>460, 1240</td><td>Wright et al., 2014</td></tr><tr><td>B-F5</td><td>460, 1240</td><td>Wright et al., 2014</td></tr><tr><td></td><td></td><td>x5</td></tr><tr><td>C</td><td>460, 560, 640, 860</td><td>Wright et al., 2014</td></tr><tr><td>C-F5</td><td>460, 560, 640, 860</td><td>Wright et al., 2014</td></tr><tr><td></td><td></td><td>x5</td></tr><tr><td>D</td><td>1240, 1640, 2120</td><td>Wright et al., 2014</td></tr><tr><td>D-F5</td><td>1240, 1640, 2120</td><td>Wright et al., 2014</td></tr><tr><td></td><td></td><td>x5</td></tr></tbody></table>
the ensemble in reproducing the evolution of an observed variable (i.e., verification variable of the truth member). The metrics obtained with and without data assimilation were compared, enabling evaluation of the impact that data assimilation had on the forecast of the ensemble. The mean of the daily Continuous Ranked Probability Score (CRPS; Hersbach, 2000; Gneiting et al., 2007; Tödter and Ahrens, 2012) over time was chosen because it is one of the most common probabilistic metrics for jointly evaluating both the reliability and sharpness of the probability distribution simulated by an ensemble system. This score decreases towards 0 when the ensemble tends towards a perfect deterministic simulation, and is expected to be reduced by data assimilation. An ensemble is reliable if events are forecast with the right probability, and has good resolution if it is able to discriminate among different events by issuing different forecasts (for further details see Atger, 1999). The reliability component of the score (RELI) can be isolated, following

Hersbach (2000). The CRPS and RELI were computed for the ensembles with and without data assimilation for the verification variables SD and SWE.

The CRPS for each time step (CRPS<sub>t</sub>) was obtained from Eq. (1)

$$CRPS_t = \int R(F_t(x) - O_t(x))^2 dx \quad (1)$$

where $F_t(x)$ is the cumulative distribution function at time t and $O_t(x)$ is the corresponding cumulative distribution function of the observation (Heviside function with center in the true value) at time t. Afterwards the CRPS<sub>t</sub> is averaged over time to obtain CRPS values of the different simulation tests. The CRPS is the sum of the reliability (RELI) of the ensemble and its resolution (RESOL)

$$CRPS = RELI + RESOL \quad (2)$$

<figure><img src="image_6.png" /><figcaption>Fig. 3. Twin experiment example of the assimilation of truth member “vii” for the 2016–17 snow season involving assimilation of the seven bands having R values from Wright et al. (2014). The upper panel shows the temporal evolution of the ensemble with and without data assimilation of the surface reflectance for the seven MODIS-like reflectances. The bottom panels show the temporal evolution of the SD (left panel graph) and the SWE (right panel graph) with and without data assimilation. P10-90 is the envelope obtained from percentiles 10 to 90 for simulations with (blue) and without (green) data assimilation. Similarly P25-75 wraps data assimilations simulations between 25 and 75 percentiles. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>
Thus, the CRPS is both a metric of ensemble reliability and its resolution. An ensemble is reliable (the lower RELI is, the better forecasts the ensemble) if events are forecasted with the right probability, and has good resolution if it is able to discriminate between different events by issuing different forecasts (Atger, 1999). For reliable systems, it is considered that the spread of the forecasted ensemble is equivalent to the resolution (Cluzet et al., 2020b), and thus it allows to easily compute the reliability.

### 3. Results

Figs. 3 and 4 show the temporal evolution of the snow surface reflectance in 2016–17 snow season for the seven bands assimilated. The open loop run is shown in green (percentiles 10–90 of this ensemble),

whereas the simulation with data assimilation is represented by blue color (percentiles 10–90 and 25–75 of the ensemble with data assimilation). As the synthetic observations of surface reflectances were not used for evaluating the simulation, and were only considered in the simulation system when a potential satellite observation was possible (see the criteria in Section 2.4), these observations are only included in the graph (red points) when the assimilation took place. Both figures include the temporal evolution of the ensemble of SD and SWE simulations, with the same color representation as for surface reflectance. Fig. 3 shows the assimilation of member “vii” (Table 1) for the 2016–17 snow season, with the seven bands assimilated using the R matrix errors of Wright et al. (2014), and Fig. 4 shows the same member, season, and bands, but with an inflation by 5 of the R matrix in the assimilation process. The dispersion of the reflectances for all bands was reduced on

<figure><img src="image_7.png" /><figcaption>Fig. 4. Twin experiment example of the assimilation of truth member “vii” for the 2016–17 snow season involving assimilation of the seven bands having R values from Wright et al. (2014) with a multiplicative factor of five. The upper panel shows the temporal evolution of the ensemble with and without data assimilation of the surface reflectance for the seven MODIS-like reflectances. The bottom panels show the temporal evolution of the SD (left panel graph) and the SWE (right panel graph) with and without data assimilation. P10-90 is the envelope obtained from percentiles 10 to 90 for simulations with (blue) and without (green) data assimilation. Similarly P25-75 wraps data assimilations simulations between 25 and 75 percentiles. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>
each assimilation date in each case. However, the assimilation involving higher R matrix errors (Fig. 4) showed slightly wider dispersion, as more members were admitted in the resampling algorithm. This was similarly observed in the dispersion of the ensemble for SD and SWE. In Fig. 3 the truth member is within the P25–75 percentile envelope from the beginning of the season until the first days of March, and within the P10–90 percentile envelope for the entire season. In Fig. 4 (higher R matrix errors) the reference values were within the P10–90 envelope, but generally not in the P25–75 envelope. This is acceptable for a data assimilation system for a relatively short time period, but should not be systematic over a longer period. Otherwise the forecast of the simulation system may advance along the different assimilations dates toward an ensemble too far from true snowpack variables, and thus diverging from ground truth.

### 3.1. Assimilation of different spectral bands combinations

The CRPS scores for both the SD (Fig. 5) and SWE (Fig. 6) for most data assimilation tests showed improved simulation results compared with the open loop simulation. A minor improvement was observed when only near infrared bands (configuration D in Fig. 5) were assimilated, as similar or slightly better CRPS scores were obtained throughout the study period for both verification variables. Of the four combinations

tested, the best result was obtained when the PF algorithm assimilated the seven spectral bands together. The use of a multiplicative factor to inflate the R matrix errors did not result in any improvement in the four band combinations tested. Despite the three snow seasons analyzed had a contrasted temporal evolution in terms of snow accumulation (Fig. 2), it was found an equivalent performance on reproducing the assimilated member with the different band combinations along the study period.

The reliability scores for our simulation system for the entire study period (Fig. 7) were lower (i.e., a more reliable ensemble) for all assimilation tests than that for the open loop ensemble. The best SD reliability value was obtained when assimilating the seven reflectance bands (configuration A), whereas for the SWE the best result was obtained for the configuration C (visible spectral bands). Given that the CRPS scores for both bulk evaluation variables improved when assimilating all spectral bands together, and that the SWE reliability improvement when assimilating visible spectral bands (configuration C) was only minor compared with that of configuration A (the seven MODIS-like reflectances), we selected this later configuration for further analysis of the assimilation of biased observations.

Using configuration A; true MODIS surface reflectance was assimilated for same dates of synthetic assimilation. The ensemble of SD obtained assimilating true satellite observations was evaluated with SD observations obtained in the automatic weather station, and thus are

<figure><img src="image_8.png" /><figcaption>**Fig. 5.** The CRPS box plots for SD for the various truth members tested for the open loop ensemble, and for the ensembles involving data assimilation simulations obtained in the various assimilation tests. These box plots include results for all spectral band combinations and the two covariance error matrixes evaluated. The top left, top right, and bottom left panels show the scores for the 2014–15, 2015–16, and 2016–17 snow seasons, respectively. The bottom right panel shows the CRPS box plot for SD for the entire study period. Black dots show the CRPS obtained in the data assimilation test of true MODIS observations, assimilating all bands (configuration with superior results). Horizontal lines inside the boxes depict mean CRPS values, boxes contain from 25 to 75 percentiles and whiskers from 10 to 90 percentiles.</figcaption></figure>
<figure><img src="image_9.png" /><figcaption>**Fig. 6.** The CRPS box plots for the SWE for the various truth members tested for the open loop ensemble, and for the ensembles involving data assimilation simulations obtained in the various assimilation tests. These box plots include results for all spectral band combinations and the two covariance error matrixes evaluated. The top left, top right, and bottom left panels show the scores for the 2014–15, 2015–16, and 2016–17 snow seasons, respectively. The bottom right panel shows the CRPS box plot for the SWE for the entire study period. Horizontal lines inside the boxes depict mean CRPS values, boxes contain from 25 to 75 percentiles and whiskers from 10 to 90 percentiles.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Fig. 7.** Reliability box plots for the entire study period with and without data assimilation for the SD (left panel) and SWE (right panel). Horizontal lines inside the boxes depict mean RELI values, boxes contain from 25 to 75 percentiles and whiskers from 10 to 90 percentiles.</figcaption></figure>

considered as ground truth. These SD represented the following percentiles of the open loop simulation: 2014–15, P68, 2015–16 P41 and 2016–17 P63. Moreover, the observed SD was nearly for the entire study period comprised in the P25-P75 envelope. The CRPS (Fig. 5) obtained for these three simulations assimilating MODIS surface reflectance

(mean CRPS 0.3 m), show that even with the bands that conveys the more information in the assimilation system, simulations results did not improve scores obtained in the open loop simulation (mean CRPS of 0.245 m). Similarly the mean reliability score for these three simulations is 0.12 m, what has same order of magnitude of that obtained for
simulations without data assimilation (Fig. 7) and nearly three times higher than the average reliability value (0.045 m) obtained assimilating synthetic reflectance with configuration A. Nonetheless, it must be highlighted that these values are obtained in three simulations along the study period and box plots in Figs. 5–7 show values obtained in 7 simulations when annual values are plotted and 21 simulations in the case of box plots for the entire study period, and thus extreme CRPS and RELI values obtained in some assimilation tests are blurred in the box-plot representation.

## 3.2. Biased surface reflectances impact

Fig. 8 shows the scores obtained when assimilating systematically biased surface reflectances during the study period. The results show that the assimilation of surface reflectances having a systematic bias of ±10% performed less well than the open loop ensemble. When the systematic bias was ±5%, similar scores to those of the open loop simulation were obtained. While our simulation system allowed a random noise of up to 5% in the assimilation reflectance, more reliable ensembles and results closer to the truth member were obtained with assimilation of reflectances involving biases <2%.

## 3.3. Impact of the reference run position within the ensemble on the assimilation efficiency

Figs. 9 and 10 show the CRPS and RELI values for the various simulations carried out for the best band configuration (all bands assimilated) and for the highest random noise in the observed reflectances that the assimilation permitted (5%). From these graphs we concluded that

the position (or percentile) of the reference run within the ensemble also impacted the score of the assimilation.

The 21 twin experiments shown in Fig. 9 demonstrated that the greatest improvement in data assimilation was obtained with smaller percentiles, while moderate improvement was obtained with medium percentiles (20–40 SD and SWE percentiles). This result does not depend on the snow season characteristics, since contrasted snow seasons in terms of snow accumulations were evaluated (Fig. 2), and this behavior (greatest improvement in data assimilation with smaller percentiles) was observed along the entire study period. Fig. 9 also shows that for some twin experiments the PF failed, and even for the best case scenario it did not improve the simulation results for the open loop ensemble. Similar conclusions were drawn for the data assimilation experiments involving a random noise between +5% and -5%.

## 4. Discussion

The assimilation of snow surface reflectance is known to improve the forecasting ability of ensemble Crocus snowpack simulations (Thirel et al., 2013; Charrois et al., 2016; Cluzet et al., 2020a), even within a spatial framework (Cluzet et al., 2020b). To build on these previous studies, our study aimed firstly to determine the optimal spectral band combination to obtain the greatest improvement in ensemble snowpack simulations. The second aim was to investigate the maximum tolerable errors in the observations acceptable for the assimilation system.

The ability to reproduce the snowpack evolution with and without data assimilation was assessed by evaluating the snowpack bulk variables SD and SWE. Overall, the assimilation of surface reflectance over three snow seasons for several dates (between 14 and 19 assimilations)

<figure><img src="image_11.png" /><figcaption>Fig. 8. CRPS (upper panels) and RELI (bottom panels) box plots for the SD (left panel) and SWE (right panel) based on assimilation of biased observations in all spectral bands. Horizontal lines inside the boxes depict mean CRPS values, boxes contain from 25 to 75 percentiles and whiskers from 10 to 90 percentiles.</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>**Fig. 9.** CRPS (upper panels) and RELI (bottom panels) versus percentile classification for the reference run for the SD (left panels) and SWE (right panels) obtained in the various twin experiments assimilating the seven MODIS-like reflectances, for simulations with (blue) and without (dark) data assimilation. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)</figcaption></figure>

improved the simulation of these two variables. The open loop ensemble of snowpack simulations was sufficiently dispersive because of the wide variety of meteorological runs and the various model parameterizations (Lafaysse et al., 2017), and this way observed SD at the automatic weather station of the study site were comprised between percentiles P41 and P68. In total, we conducted 276 assimilation tests based on the random extraction of 21 truth members from the open loop ensemble. We evaluated the impact of the assimilation of members having intermediate, deep, and shallow snowpacks. The percentile classification of the truth members showed that greater improvement was obtained for shallow SD (Figs. 9 and 10), with satisfactory improvement of the ensemble simulations for most of the assimilation tests. Probably the simulation of shallow snow accumulations is more sensible to small differences when the PF selects (and replicates) those members that are closer to the observed surface reflectance, than when simulating a thick snowpack. Oppositely, the three ensemble simulations assimilating true MODIS surface reflectance observations did not improve simulation results.

The number of simulations that comprise the ensemble is an important issue when simulating over large domains, because an exponential number of in situ simulations is required (as shown in the example of the Col du Lautaret), which markedly increases the requirement for computational resources. Larger ensembles have been shown to be associated with a major decrease in forecast errors when simulating the snowpack (e.g., 2000 members; Magnusson et al., 2017). Similarly, smaller ensembles including 300 members (Charrois et al.,

2016) and 100 members (Piazzi et al., 2018) have also shown good results and improved simulation performance. The latter study showed that taking into account the uncertainty inherent in the model using a stochastic procedure ensures a suitable spread of the ensemble with a moderate number of particles (100 in their case). Consequently, we followed a similar procedure to that of Piazzi et al. (2018) and proposed by (Moradkhani et al., 2005), but took advantage of the Crocus multi-physics ensemble ESCROC (Lafaysse et al., 2017). For each assimilation step, following resampling of particles through the SIR algorithm both the ensemble of meteorological forcings and the various ESCROC configurations were randomly assigned to the particles issued by the filter. We showed that this procedure was effective with our ensemble of 200 members for simulating snowpack evolution.

With respect to the various combinations of spectral bands, the CRPS for both the SD and SWE for most combinations was smaller than that obtained for the open loop ensemble. For the entire study period, only configuration D (near-infrared reflectance only) showed higher CRPS scores than those of the reference run versus the open loop ensemble.

The greater improvement obtained when assimilating the first seven MODIS-like bands indicates that information over the entire spectrum is required in the assimilation system. Visible bands are known to provide information about the LAP content of the snow surface, whereas near-infrared bands are more sensitive to snow grain SSA. During the snow season, both LAP deposition and changes in the SSA occur. Therefore, assimilating information on the snow surface reflectance in various regions of the solar spectrum is strongly recommended.
<figure><img src="image_13.png" /><figcaption>**Fig. 10.** CRPS (upper panels) and RELI (bottom panels) versus percentile classification for the reference run for the SD (left panels) and SWE (right panels) obtained in the various twin experiments assimilating the seven MODIS-like reflectances with a random noise of ±5%, for simulations with (blue) and without (dark) data assimilation.</figcaption></figure>

However, some data assimilation twin experiments failed, including for the best case configuration (assimilation of all MODIS-like reflectance bands), with similar or worse results than for the open loop simulation. These failures in the assimilation system have two plausible explanations. One is degeneration of the ensemble to too few particles, which requires the use of strategies to tackle and mitigate this problem (Cluzet et al., 2020b). Another explanation is the nature of the observations assimilated. For instance, when assimilating true MODIS surface reflectance, if after a long period without solid precipitation (e.g., during the melt period) a small amount of snowfall accumulates some hours prior to a valid satellite acquisition, and is not adequately simulated by the meteorological forcing, the snow surface reflectance describing the snowpack will differ substantially from the simulated value. Similarly, an erroneous simulation of surface reflectance when assimilating synthetic observations, may introduce deviation which potentially originate an ensemble far from the reference run. Under these circumstances, the PF algorithm would select particles not accurately simulated, contributing to degeneration of the ensemble.

Our PF data assimilation scheme failed to improve simulation results when assimilating true MODIS snow surface reflectance. This may be related to the errors in surface reflectance retrieved from MODIS in complex terrain. The study site is locally flat close to the meteorological station but the surface reflectance retrieved by the space borne MODIS sensor is affected by complex terrain effects (e.g. Lamare et al., 2020). The small deviations between ground based and satellite observations

reported in Greenland by Wright et al., (2014) are likely not transferable to such sites in complex terrain. Snow observations from optical satellite sensors are indeed affected by various shortcomings, mainly because of the complex interactions between radiation and topography (Dumont and Gascoin 2016; Lamare et al., 2020); these often result in biased observations (Cluzet et al., 2020a). However, to the best of our knowledge, the accuracy of MODIS spectral surface reflectance in highly heterogeneous mountain areas has not been accurately determined, and thus accuracy assessments from other areas having distinct characteristics are applied. Several studies have attempted to use MODIS observations in complex topography but the assimilation experiment were either non successful (Charrois et al., 2016) or non-possible because of the high observation errors (Cluzet et al., 2020a). On the contrary, the assimilation of ground based snow surface reflectance observations using the same spectral bands as MODIS improved simulation results (Charrois 2017). Therefore, assimilating accurate snow surface reflectance guarantees more reliable snowpack simulations.

Aiming to provide guidelines about the maximum deviation from ground truth that the assimilation of snow surface reflectance admits to improve simulations, differing biases to the synthetic observations were tested. We showed that the maximum deviation from ground truth that our data assimilation scheme could tolerate before there was a decrease in the simulation performance was ±5% for white noise and ±2% for systematic bias. Although these criteria could be met for flat snow-covered surfaces (Kokhanovsky et al., 2019; Moustafa et al., 2017),
they are beyond the current limit of standard MODIS products in heterogeneous mountain areas (Masson et al., 2018; Cluzet et al., 2020a). However, recent progress in accounting for topographic effects (Lamare et al., 2020), combined with increased availability of higher resolution satellite data (Gascoin et al., 2019), suggest that such criteria may soon be met, even in mountainous areas. Our findings confirm that when remote sensing science can achieve the deviations from ground truth reported here, the assimilation of data from optical satellite sensors will improve ensemble forecasting capabilities.

## 5. Conclusions

Our results demonstrate that the assimilation of snowpack surface reflectances using the Particle Filter algorithm improves simulation of the temporal evolution of snowpack bulk variables. Our study used synthetic observations, but the assimilation dates were selected from true MODIS images made under appropriate conditions (clear sky, full snow cover, viewing zenith angle >60°) for retrieving surface reflectance from the snowpack. These restricted conditions for acquiring snow surface reflectance led to the assimilation of data from 14 to 20 dates in each of three snow seasons. Despite the reduced number of assimilation dates over more than seven months of snowpack simulations, the simulation system showed an improvement in results compared with the simulations made without data assimilation. The assimilation of the first seven MODIS-like bands, including information in visible and near-infrared wavelengths, showed the best results from all band combinations tested. In the contrary, the assimilation of true MODIS surface reflectance failed and did not improve open loop simulations results. In light of this shortcoming, different biases were introduced to the synthetic observations to obtain the maximum deviation from ground truth that data assimilation admits. The simulation system enabled the assimilation of surface reflectances having a random white noise less than or equal to 5%, or a systematic bias less than or equal to ±2% from ground truth. This indicates the maximum deviation that satellite-derived products will have to provide for improvement of ensemble simulations through data assimilation to occur. The results of this study provide the specifications in terms of accuracy for the processing of satellite reflectances.

## CRediT authorship contribution statement

**J. Revuelto:** Methodology, Data curation, Formal analysis, Resources, Investigation, Validation, Visualization, Writing - original draft, Writing - review & editing, Project administration, Funding acquisition.
**B. Cluzet:** Methodology, Software, Validation, Data curation, Investigation, Writing - review & editing. **N. Duran:** Software, Validation, Visualization. **M. Fructus:** Software, Validation, Visualization. **M. Lafaysse:** Conceptualization, Software, Writing - review & editing. **E. Cosme:** Conceptualization, Writing - review & editing. **M. Dumont:** Conceptualization, Methodology, Supervision, Writing - review & editing, Project administration, Funding acquisition.

## Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Acknowledgments

J. Revuelto was supported by a Post-doctoral Fellowship of the AXA research fund (le Post-Doctorant Jesús Revuelto est bénéficiaire d'une bourse postdoctorale du Fonds AXA pour la Recherchem Ref: CNRM 3.2.01/17). J. Revuelto is now supported by a "Juan de la Cierva Incorporación" Fellowship of the Spanish Ministry of Science and Innovation (IJC2018-036260-I). IGE and CNRM/CEN are part of Labex

OSUG@2020. This work was partly supported by the French national program LEFE/INSU (ASSURANCE), APR CNES MIOSOTI, and an ANR JCJC EBONI grant (ANR-16-CE01-0006). This research was also supported by Lautaret Garden-UMS 3370 (Univ. Grenoble Alpes, CNRS, SAJF, 38000 Grenoble, France), a member of AnaEE-France (ANR-11-INBS-0001AnaEE-Services, Investissements d'Avenir frame) and the eLTER-Europe network (Univ. Grenoble Alpes, CNRS, LSTER Zone Atelier Alpes, 38000 Grenoble, France).

## References

* Andreadis, K.M., Lettenmaier, D.P., 2006. Assimilating remotely sensed snow observations into a macroscale hydrology model. Adv. Water Resour. 29 (6), 872-886.

* Atger, F., 1999. The skill of ensemble prediction systems. Mon. Weather Rev. 127 (9), 1941-1953.

* Baba, M., Gascoin, S., Hanich, L., 2018. Assimilation of Sentinel-2 Data into a Snowpack Model in the High Atlas of Morocco. Remote Sensing 10 (12), 1982. https://doi.org/10.3390/rs10121982.

* Browne, P.A., van Leeuwen, P.J., 2015. Twin experiments with the equivalent weights particle filter and HadCM3. Q. J. R. Meteorol. Soc. 141 (693), 3399-3414.

* Charrois, L. (2017). PhD dissertation : Assimilation de réflectances satellitaires du domaine visible en proche infrarouge dans un modèle détaillé de manteau neigeux. Sciences de la Terre. Université Grenoble Alpes, 2017. Français. (NNT : 2017GREAU001). (tel-01492360v2).

* Charrois, L., Cosme, E., Dumont, M., Lafaysse, M., Morin, S., Libois, Q., Picard, G., 2016. On the assimilation of optical reflectances and snow depth observations into a detailed snowpack model. Cryosphere 10 (3), 1021-1038.

* Clark, M.P., Rupp, D.E., Woods, R.A., Zheng, X., Ibbitt, R.P., Slater, A.G., Schmidt, J., Uddstrom, M.J., 2008. Hydrological data assimilation with the ensemble Kalman filter: Use of streamflow observations to update states in a distributed hydrological model. Adv. Water Resour. 31 (10), 1309-1324.

* Cluzet, B., Revuelto, J., Lafaysse, M., Tuzet, F., Cosme, E., Picard, G., Arnaud, L., Dumont, M., 2020a. Towards the assimilation of satellite reflectance into semi-distributed ensemble snowpack simulations. Cold Reg. Sci. Technol. 170, 102918. https://doi.org/10.1016/j.coldregions.2019.102918.

* Cluzet, B., Lafaysse, M., Cosme, E., Albergel, C., Meunier, L.-F., Dumont, M., 2020b. CrocO v1.0: a Particle Filter to assimilate snowpack observations in a spatialised framework. Geoscientific Model Development Discussions 1-36.

* Cortés, G., Girotto, M., Margulis, S., 2016. Snow process estimation over the extratropical Andes using a data assimilation framework integrating MERRA data and Landsat imagery. Water Resour. Res. 52 (4), 2582-2600.

* Decharme, B., Brun, E., Boone, A., Delire, C., Le Moigne, P., Morin, S., 2016. Impacts of snow and organic soils parameterization on northern Eurasian soil temperature profiles simulated by the ISBA land surface model. Cryosphere 10 (2), 853-877.

* Dozier, J., Green, R.O., Nolin, A.W., Painter, T.H., 2009. Interpretation of snow properties from imaging spectrometry. Remote Sens. Environ. 113, S25-S37.

* Dubinkina, S., Goosse, H., Sallaz-Damaz, Y., Crespin, E., and Crucifix, M. (2011). Testing a particle filter to reconstruct climate changes over the past centuries. Int. J. Bifurc. Chaos 21, 3611-3618.

* Dumont, M., Durand, Y., Arnaud, Y., Six, D., 2012. Variational assimilation of albedo in a snowpack model and reconstruction of the spatial mass-balance distribution of an alpine glacier. J. Glaciology 58, 151-164.

* Dumont, M., Gascoin, S., 2016. 4 - Optical Remote Remote Sensing of Snow Cover. In Land Surface Remote Sensing in Continental Hydrology, N. Baghdadi, and M. Zribi, eds. (Elsevier), pp. 115-137.

* Durand, Y., Latenser, M., Giraud, G., Etchevers, P., Lesaffre, B., Mérdinol, L., 2009. Reanalysis of 44 yr of climate in the French Alps (1958-2002): methodology, model validation, climatology, and trends for air temperature and precipitation. J. Appl. Meteor. Climatol. 48, 429-449.

* Evensen, G., 2003. The Ensemble Kalman Filter: theoretical formulation and practical implementation. Ocean Dyn. 53 (4), 343-367.

* Frei, A., Tedesco, M., Lee, S., Foster, J., Hall, D.K., Kelly, R., Robinson, D.A., 2012. A review of global satellite-derived snow products. Adv. Space Res. 50 (8), 1007-1029.

* Gaál, L., Szolgay, J., Kohnová, S., Hlavcová, K., Parajka, J., Viglione, A., Merz, R., Blöschl, G., 2015. Dependence between flood peaks and volumes: a case study on climate and hydrological controls. Hydrol. Sci. J. 60 (6), 968-984.

* Gascoin, S., Hagolle, O., Huc, M., Jarlan, L., Dejoux, J.-F., Szczypta, C., Marti, R., Sánchez, R., 2015. A snow cover climatology for the Pyrenees from MODIS snow products. Hydrol. Earth Syst. Sci. 19 (5), 2337-2351.

* Gascoin, S., Grizonnet, M., Bouchet, M., Salgues, G., Hagolle, O., 2019. Theia Snow collection: high-resolution operational snow cover maps from Sentinel-2 and Landsat-8 data. Earth Syst. Sci. Data 11 (2), 493-514.

* Girotto, M., Musselman, K.N., Essery, R.L.H., 2020. Data assimilation improves estimates of climate-sensitive seasonal snow. Curr. Clim. Change Rep. 6 (3), 81-94.

* Gneiting, T., Balabdaoui, F., Raftery, A.E., 2007. Probabilistic forecasts, calibration and sharpness. J. R. Stat. Soc. Ser. B Stat. Methodol. 69 (2), 243-268.

* Gordon, N.J., Salmon, D.J., Smith, A.F.M., 1993. Novel approach to nonlinear/non-Gaussian Bayesian state estimation. IEE Proc. F Radar Signal Process. 140 (2), 107. https://doi.org/10.1049/ip-f-2.1993.0015.

* Hall, D.K., Riggs, G.A., 2007. Accuracy assessment of the MODIS snow products. Hydrol. Process. 21 (12), 1534-1547.
Hersbach, H., 2000. Decomposition of the continuous ranked probability score for ensemble prediction systems. Weather Forecast. 15 (5), 559-570.

Immezerle, W.W., Droogers, P., de Jong, S.M., Bierkens, M.F.P., 2009. Large-scale monitoring of snow cover and runoff simulation in Himalayan river basins using remote sensing. Remote Sens. Environ. 113 (1), 40-49.

Janjić, T., Bormann, N., Bocquet, M., Carton, J.A., Cohn, S.E., Dance, S.L., Losa, S.N., Nichols, N.K., Potthast, R., Waller, J.A., Weston, P., 2018. On the representation error in data assimilation. Quart. J. R. Meteorol. Soc. 144 (713), 1257-1278.

Josse, B., Simon, P., Peuch, V.-H., 2004. Radon global simulations with the multiscale chemistry and transport model MOCAGE. Tellus B Chem. Phys. Meteorol. 56 (4), 339-356.

Kokhanovsky, A., Lamare, M., Danne, O., Brockmann, C., Dumont, M., Picard, G., Arnaud, L., Favier, V., Jourdain, B., Le Meur, E., Di Mauro, B., Aoki, T., Niwano, M., Rozanov, V., Korkin, S., Kipfstuhl, S., Freitag, J., Hoerhold, M., Zuhr, A., Vladimirova, D., Faber, A.-K., Steen-Larsen, H., Wahl, S., Andersen, J., Vandecrux, B., van As, D., Mankoff, K., Kern, M., Zege, E., Box, J., 2019. Retrieval of snow properties from the Sentinel-3 ocean and land colour instrument. Remote Sensing 11 (19), 2280. https://doi.org/10.3390/rs11192280.

Lafaysse, M., Cluzet, B., Dumont, M., Lejeune, Y., Vionnet, V., Morin, S., 2017. A multispatial ensemble system of numerical snow modelling. Cryosphere 11 (3), 1173-1198.

Lamare, M., Dumont, M., Picard, G., Larue, F., Tuzet, F., Delcourt, C., Arnaud, L., 2020. Simulating optical top-of-atmosphere radiance satellite images over snow-covered rugged terrain. Cryosphere 14 (11), 3995-4020.

Larue, F., Royer, A., De Sève, D., Roy, A., Cosme, E., 2018. Assimilation of passive microwave AMSR-2 satellite observations in a snowpack evolution model over northeastern Canada. Hydrol. Earth Syst. Sci. 22 (11), 5711-5734.

Larue, F., Picard, G., Arnaud, L., Ollivier, I., Delcourt, C., Lamare, M., Tuzet, F., Revuelto, J., Dumont, M., 2020. Snow albedo sensitivity to macroscopic surface roughness using a new ray-tracing model. Cryosphere 14 (5), 1651-1672.

Largeron, C., Dumont, M., Morin, S., Boone, A., Lafaysse, M., Metref, S., Cosme, E., Jonas, T., Winstral, M., Margulis, S.A., 2020. Toward snow cover estimation in mountainous areas using modern data assimilation methods: a review. Front. Earth Sci. 8.

van Leeuwen, P.J., 2009. Particle filtering in geophysical systems. Mon. Weather Rev. 137, 4089-4114.

Lehning, M., Völksch, I., Gustafsson, D., Nguyen, T.A., Stähli, M., Zappa, M., 2006. ALPINE3D: a detailed model of mountain surface processes and its application to snow hydrology. Hydrol. Process. 20 (10), 2111-2128.

Libois, Q., Picard, G., Arnaud, L., Dumont, M., Lafaysse, M., Morin, S., Lefebvre, E., 2015. Summertime evolution of snow specific surface area close to the surface on the Antarctic Plateau. The Cryosphere 9 (6), 2383-2398.

Magnusson, J., Winstral, A., Stordal, A.S., Essery, R., Jonas, T., 2017. Improving physically based snow simulations by assimilating snow depths using the particle filter. Water Resour. Res. 53 (2), 1125-1143.

Margulis, S.A., Fang, Y., Li, D., Lettenmaier, D.P., Andreadis, K., 2019. The utility of infrequent snow depth images for deriving continuous space-time estimates of seasonal snow water equivalent. Geophys. Res. Lett. 46 (10), 5331-5340.

Masson, T., Dumont, M., Mura, M., Sirguey, P., Gascoin, S., Dedieu, J.-P., Chanussot, J., 2018. An assessment of existing methodologies to retrieve snow cover fraction from MODIS data. Remote Sensing 10 (4), 619. https://doi.org/10.3390/rs10040619.

Matgen, P., Montanari, M., Hostache, R., Pfister, L., Hoffmann, L., Plaza, D., Pauwels, V. R.N., De Lannoy, G.J.M., De Keyser, R., Savenije, H.H.G., 2010. Towards the sequential assimilation of SAR-derived water stages into hydraulic models using the Particle Filter: proof of concept. Hydrol. Earth Syst. Sci. 14 (9), 1773-1785.

Ménard, C.B., Essery, R., Barr, A., Bartlett, P., Derry, J., Dumont, M., Fierz, C., Kim, H., Kontu, A., Lejeune, Y., Marks, D., Niwano, M., Raleigh, M., Wang, L., Wever, N., 2019. Meteorological and evaluation datasets for snow modelling at 10 reference sites: description of in situ and bias-corrected reanalysis data. Earth Syst. Sci. Data 11 (2), 865-880.

Moradkhani, H., Hsu, K.-L., Gupta, H., Sorooshian, S., 2005. Uncertainty assessment of hydrologic model states and parameters: Sequential data assimilation using the particle filter. Water Resour. Res. 41 (5) https://doi.org/10.1029/2004WR003604.

Morin, S., Horton, S., Techel, F., Bavay, M., Coléou, C., Fierz, C., Gobiet, A., Hagenmuller, P., Lafaysse, M., Lizar, M., Mitterer, C., Monti, F., Müller, K., Olefs, M., Snook, J.S., van Herwijnen, A., Vionnet, V., 2020. Application of physical snowpack models in support of operational avalanche hazard forecasting: A status report on current implementations and prospects for the future. Cold Reg. Sci. Technol. 170, 102910. https://doi.org/10.1016/j.coldregions.2019.102910.

Moustafa, S.E., Rennermalm, A.K., Roman, M.O., Wang, Z., Schaaf, C.B., Smith, L.C., Koenig, L.S., Erb, A., 2017. Evaluation of satellite remote sensing albedo retrievals over the ablation area of the southwestern Greenland ice sheet. Remote Sens. Environ. 198, 115-125.

Painter, T.H., Bryant, A.C., Skiles, S.M., 2012. Radiative forcing by light absorbing impurities in snow from MODIS surface reflectance data. Geophys. Res. Lett. 39 (17), n/a-n/a.

Parajka, J., Blöschl, G., 2008. The value of MODIS snow cover data in validating and calibrating conceptual hydrologic models. Journal of Hydrology 358, 240-258.

Piazzì, G., Thirel, G., Campo, L., Gabellani, S., 2018. A particle filter scheme for multivariate data assimilation into a point-scale snowpack model in an Alpine environment. Cryosphere 12 (7), 2287-2306.

Piazzì, G., Campo, L., Gabellani, S., Castelli, F., Cremonese, E., Cella, U.M. di Stevenin, H., Ratto, S.M., 2019. An Enkf-Based Scheme for Snow Multivariable Data Assimilation at an Alpine Site. J. Hydrol. Hydromech. 67, 4-19.

Revuelto, J., Lecourt, M., Lafaysse, M., Zin, I., Charrois, L., Vionnet, V., Dumont, M., Rabatel, A., Six, D., Condom, T., Morin, S., Viani, A., Sirguey, P., 2018. Multi-criteria evaluation of snowpack simulations in complex alpine terrain using satellite and in situ observations. Remote Sens. 10 (8), 1171. https://doi.org/10.3390/rs10081171.

Revuelto, J., Billecocq, P., Tuzet, F., Cluzet, B., Lamare, M., Larue, F., Dumont, M., 2020. Random forests as a tool to understand the snow depth distribution and its evolution in mountain areas. Hydrol. Process. 34 (26), 5384-5401.

Schweizer, J., Jamieson, J.B., Schneebeli, M., 2003. Snow avalanche formation. Rev. Geophys. 41, 2.1-2.25.

Schweizer, J., Kronholm, K., Jamieson, J.B., Birkeland, K.W., 2008. Review of spatial variability of snowpack properties and its importance for avalanche formation. Cold Reg. Sci. Technol. 51 (2-3), 253-272.

Scipión, D.E., Mott, R., Lehning, M., Schneebeli, M., Berne, A., 2013. Seasonal small-scale spatial variability in alpine snowfall and snow accumulation. Water Resour. Res. 49 (3), 1446-1457.

Seidel, F.C., Rittger, K., Skiles, S.M., Molotch, N.P., Painter, T.H., 2016. Case study of spatial and temporal variability of snow cover, grain size, albedo and radiative forcing in the Sierra Nevada and Rocky Mountain snowpack derived from imaging spectroscopy. Cryosphere 10 (3), 1229-1244.

Sirguey, P., Mathieu, R., Arnaud, Y., 2009. Subpixel monitoring of the seasonal snow cover with MODIS at 250 m spatial resolution in the Southern Alps of New Zealand: Methodology and accuracy assessment. Remote Sens. Environ. 113, 160-181.

Sirguey, P., Mathieu, R., Arnaud, Y., Khan, M.M., Chanussot, J., 2008. Improving MODIS spatial resolution for snow mapping using wavelet fusion and ARSIS concept. IEEE Geoscience and Remote Sensing Letters 5, 78-82.

Sirguey, P., Still, H., Cullen, N.J., Dumont, M., Arnaud, Y., Conway, J.P., 2016. Reconstructing the mass balance of Brewster Glacier, New Zealand, using MODIS-derived glacier-wide albedo. Cryosphere 10 (5), 2465-2484.

Skiles, S.M., Flanner, M., Cook, J.M., Dumont, M., Painter, T.H., 2018. Radiative forcing by light-absorbing particles in snow. Nat. Clim. Change 1.

Smyth, E.J., Raleigh, M.S., Small, E.E., 2019. Particle filter data assimilation of monthly snow depth observations improves estimation of snow density and SWE. Water Resour. Res. 55 (2), 1296-1311.

Swinbank, R., Kyouda, M., Buchanan, P., Froude, L., Hamill, T.M., Hewson, T.D., Keller, J.H., Matsueda, M., Methven, J., Pappenberger, F., et al., 2016. The TIGGE Project and Its Achievements. Bull. Amer. Meteor. Soc. 97, 49-67.

Thirel, G., Salamon, P., Burek, P., Kalas, M., 2013. Assimilation of MODIS snow cover area data in a distributed hydrological model using the particle filter. Remote Sens. 5 (11), 5825-5850.

Tödter, J., Ahrens, B., 2012. Generalization of the ignorance score: continuous ranked version and its decomposition. Mon. Weather Rev. 140, 2005-2017.

Tuzet, F., Dumont, M., Lafaysse, M., Picard, G., Arnaud, L., Voisin, D., Lejeune, Y., Charrois, L., Nabat, P., Morin, S., 2017. A multilayer physically based snowpack model simulating direct and indirect radiative impacts of light-absorbing impurities in snow. Cryosphere 11 (6), 2633-2653.

Tuzet, F., Dumont, M., Picard, G., Lamare, M., Voisin, D., Nabat, P., Lafaysse, M., Larue, F., Revuelto, J., Arnaud, L., 2020. Quantification of the radiative impact of light-absorbing particles during two contrasted snow seasons at Col du Lautaret (2058 m.a.s.l., French Alps). Cryosphere Discussions 1-38.

Vernay, M., Lafaysse, M., Mérdindol, L., Giraud, G., Morin, S., 2015. Ensemble forecasting of snowpack conditions and avalanche hazard. Cold Reg. Sci. Technol. C 120, 251-262.

Vionnet, V., Brun, E., Morin, S., Boone, A., Faroux, S., Le Moigne, P., Martin, E., Willemet, J.-M., 2012. The detailed snowpack scheme Crocus and its implementation in SURFEX v7.2. Geosci. Model Dev. 5 (3), 773-791.

Warren, S.G., 1982. Optical properties of snow. Rev. Geophys. 20 (1), 67. https://doi.org/10.1029/RG020i001p00067.

Wayand, N.E., Lundquist, J.D., Clark, M.P., 2015. Modeling the influence of hypsometry, vegetation, and storm energy on snowmelt contributions to basins during rain-on-snow floods. Water Resour. Res. 51 (10), 8551-8569.

Winstral, A., Magnusson, J., Schirmer, M., Jonas, T., 2019. The bias-detecting ensemble: a new and efficient technique for dynamically incorporating observations into physics-based. Multilayer Snow Models. Water Resour. Res. 55 (1), 613-631.

Wright, P., Bergin, M., Dibb, J., Lefer, B., Domine, F., Carman, T., Carmagnola, C., Dumont, M., Courville, Z., Schaaf, C., Wang, Z., 2014. Comparing MODIS daily snow albedo to spectral albedo field measurements in Central Greenland. Remote Sens. Environ. 140, 118-129.