# Analyzing the sensitivity of a blowing snow model (SnowPappus) to precipitation forcing, blowing snow, and spatial resolution

Ange Haddjeri, Matthieu Baron, Matthieu Lafaysse, Louis Le Toumelin, César Deschamps-Berger, Vincent Vionnet, Simon Gascoin, Matthieu Vernay, Marie Dumont

► To cite this version:

Ange Haddjeri, Matthieu Baron, Matthieu Lafaysse, Louis Le Toumelin, César Deschamps-Berger, et al. Analyzing the sensitivity of a blowing snow model (SnowPappus) to precipitation forcing, blowing snow, and spatial resolution. The Cryosphere, 2024, 18 (7), pp.3081-3116. 10.5194/tc-18-3081-2024. hal-04634649

HAL Id: hal-04634649

https://hal.science/hal-04634649v1

Submitted on 4 Jul 2024

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
# Analyzing the sensitivity of a blowing snow model (SnowPappus) to precipitation forcing, blowing snow, and spatial resolution

Ange Haddjeri¹, Matthieu Baron¹, Matthieu Lafaysse¹, Louis Le Toumelin¹, César Deschamps-Berger²,³, Vincent Vionnet⁴, Simon Gascoin³, Matthieu Vernay¹, and Marie Dumont¹

¹Univ. Grenoble Alpes, Université de Toulouse, Météo-France, CNRS, CNRM, Centre d'Études de la Neige, Grenoble, France

²Instituto Pirenaico de Ecología, Consejo Superior de Investigaciones Científicas (IPE-CSIC), Zaragoza, Spain

³Centre d'Études Spatiales de la Biosphère, CESBIO, Univ. Toulouse, CNES/CNRS/INRAE/IRD/UPS, Toulouse, France

⁴Meteorological Research Division, Environment and Climate Change Canada, Dorval, QC, Canada

**Correspondence:** Ange Haddjeri (ange.haddjeri@umr-cnrm.fr)

Received: 3 November 2023 – Discussion started: 23 November 2023

Revised: 3 May 2024 – Accepted: 10 May 2024 – Published: 4 July 2024

**Abstract.** Accurate snow cover modeling is a high-stakes issue for mountain regions. Alpine snow evolution and spatial variability result from a multitude of complex processes including interactions between wind and snow. The SnowPappus blowing snow model was designed to add blowing snow modeling capabilities to the SURFEX/Crocus simulation system for applications across large spatial and temporal extents. This paper presents the first spatialized evaluation of this simulation system over a 902 km² domain in the French Alps. Here we compare snow cover simulations to the spatial distribution of snow height obtained from Pléiades satellites stereo imagery and to snow melt-out dates from Sentinel-2 and Landsat 8 time series over three snow seasons. We analyzed the sensitivity of the simulations to three different precipitation datasets and two horizontal resolutions. The evaluations are presented as a function of elevation and landform types. The results show that the SnowPappus model forced with high-resolution wind fields enhances the snow cover spatial variability at high elevations allowing a better agreement between observed and simulated spatial distributions above 2500 m and near peaks and ridges. Model improvements are not obvious at low to medium altitudes where precipitation errors are the prevailing uncertainty. Our study illustrates the necessity to consider error contributions from blowing snow, precipitation forcings, and unresolved subgrid variability for robust evaluations of spatialized snow simulations. Despite the significant effect of the unresolved spatial scales of snow transport, 250 m horizontal-resolution snow

simulations using SnowPappus are found to be a promising avenue for large-scale modeling of alpine snowpacks.

## 1 Introduction

Snow cover in mountainous terrains is characterized by an important variability at multiple temporal and spatial scales (Pomeroy and Gray, 1995; Clark et al., 2011; Anderson et al., 2014; Mott et al., 2018). This variability results from a large diversity of slopes, aspects, and elevations which induce a high variability in precipitation, wind, temperature, and radiation. Accurate snowpack modeling is key to describing this high variability in hydrological applications, climate projections, and hazard forecasting in mountainous terrains (IPCC, 2022; Morin et al., 2020).

The spatial variability in precipitation amount and phase at different scales is known to be one of the main sources of snowpack variability at the mountain range scale (1–100 km) (Clark et al., 2011). These patterns are specific to the regional topography and atmospheric flow, with regions of increased or decreased precipitation (Colle et al., 2013). At the slope scale (a few hundred meters), preferential deposition and snowfall enhancement are the main processes responsible for snowfall variability (before snowflake settlement) (Mott et al., 2018). The preferential deposition is the result of the interaction of the near-surface flow field with particle trajectories, creating areas of local accumulation (Lehning
et al., 2008). Snowfall enhancement is a process in which surface flows are responsible for locally increased air moisture, leading to the formation or maintenance of low-level clouds, ultimately enhancing solid precipitations through a seeder-feeder mechanism (Bergeron, 1965; Choularton and Perry, 1986; Minder et al., 2011). Finally, post-depositional processes determine the snow variability at the slope scale and below. Blowing snow transport, snow sublimation, snow redistribution by avalanches, snow compaction, and melt are the main processes at work (Winstral et al., 2002; Bernhardt and Schulz, 2010; Mott et al., 2018). Blowing snow transport produces a mass transfer of snow from windward areas to leeward deposition zones. Additional mass loss occurs in those events due to the sublimation of suspended snow (Liston and Sturm, 1998; Yang et al., 2010). The heterogeneity of snow mass loss due to snowmelt (Brauchli et al., 2017) and surface sublimation (Pomeroy et al., 1998; Strasser et al., 2008) also contributes to snowpack variability in alpine terrain.

The multi-scale variability in the alpine snowpack makes it challenging to represent in numerical models. Various approaches have been developed for that purpose. Numerical weather prediction (NWP) models at the kilometer scale or higher-resolution models can explicitly represent the orographic precipitation patterns and in some cases local low-level cloud formation, triggering snowfall enhancement (Lehning et al., 2006; Vionnet et al., 2017; Wang and Huang, 2017; Monteiro et al., 2022). However, the quantification of the precipitation is still impacted by important uncertainties with errors in precipitation amounts and phases, localization, and timing increased with coarse grid size (Clark et al., 2011; Ménard et al., 2019; Lundquist et al., 2019). Downscaling tools can be used to better represent the local meteorology from NWP models (Sen Gupta and Tarboton, 2016; Marsh et al., 2023; Bernhardt et al., 2010; Mital et al., 2022).

Due to the complex intertwining of snow variability processes, in particular blowing snow transport, the community found benefits in the development of dedicated high-resolution models coupled to (e.g., Vionnet et al., 2014; Sharma et al., 2023) or forced by atmospheric models (e.g., Lehning et al., 2006; Liston et al., 2007; Marsh et al., 2020; Baron et al., 2024; Quéno et al., 2023). The spatial evaluation of this type of system, dedicated to modeling part of the observed snow spatial variability, is a challenge in itself. The snow modeling community has long been evaluating this kind of model locally. Evaluations can be carried out using direct measurements of the variable of interest; for example for blowing snow modeling, Vionnet et al. (2014) evaluated directly the simulated blowing snow fluxes of the Meso-NH/Crocus system against locally measured fluxes from a snow particle counter (Sato et al., 1993). Amory et al. (2021) also compare simulated snow transport occurrence and mass fluxes to on-field observations in Antarctica. Similarly, Baron et al. (2024) compare the simulated snow transport occurrence and transport mass fluxes from the SnowPapp-

pus model to field observations. One of the main drawbacks of this method is the low number of direct snow transport observations found in the literature, combined with the very high spatial variability in snow transport fluxes. As a result, spatialized snow simulations are more classically evaluated using measurements of snow height and snow water equivalent (SWE). Prasad et al. (2001) and Liston et al. (2007) compare snow simulations done using the SnowTran-3D system (Liston and Sturm, 1998; Liston et al., 2007) with SWE measurements from manual snow surveys. Marsh et al. (2020) compare simulated SWE from the Canadian Hydrological Model (CHM) to point observations. Mott et al. (2008) evaluate the Alpine3D system (Lehning et al., 2006), comparing simulation with an interpolated map of snow height measurements. The use of interpolation methods in Mott et al. (2008) is justified by the high measurement density of the study area. However, in the majority of other study areas, in situ measurements are too sparse to characterize the complex spatial variability in alpine snow cover (Pepin et al., 2015; Bales et al., 2006; Vernay et al., 2022; Pomeroy et al., 2009). This situation explains the recent use of snow remote sensing (uncrewed aerial vehicle, satellite) as a means of spatial evaluation for large-scale snow cover models. Those methods are still limited by the low availability of remotely sensed variables, mainly snow cover fraction and derived variables or snow height. For example, Vionnet et al. (2021) compare CHM snow simulations (Marsh et al., 2020) with airborne Lidar snow height maps as well as Sentinel-2 snow cover maps. Very recently, Quéno et al. (2023) used Lidar-acquired snow height maps to evaluate the ability of the FSM2oshd framework (Quéno et al., 2023) to represent snow accumulation and erosion areas.

In the above-cited studies using the snow height, SWE, or snow presence variables, the evaluation experiments are generally carried out using a single set of meteorological inputs. However, it is known that meteorological inputs explain an important part of the variability (Clark et al., 2011; Colle et al., 2013; Mott et al., 2018) and uncertainty (Raleigh et al., 2015; Günther et al., 2019) in simulated snow height and SWE. While the evaluations of model inputs and simulated processes are usually conducted separately, they are often interdependent. Thus, it is difficult to determine if potential errors come from the model or the input variable. For instance, a good simulation of snow height requires both accurate precipitation input and a robust snow evolution model. To robustly assess the value of a distributed snow model, the uncertainty in meteorological forcing must be considered.

To improve the snow spatial variability in the French snow modeling system, Baron et al. (2024) have developed a novel, explicit, blowing snow transport model, SnowPappus, coupled with the Crocus physical snow simulation model (Vionnet et al., 2012). This model explicitly represents the vertically integrated saltation and suspension mass fluxes as a function of wind speed and surface snowpack properties simulated by Crocus. This system is implemented on a
regular grid to simulate snowpack evolution at a mountain range scale (approx. 100 000 km²), for multiple snow seasons. Baron et al. (2024) focused on the evaluation of the model at point scale, demonstrating its ability to accurately simulate blowing snow fluxes and occurrence at observation stations. However, spatial evaluations of simulations on large distributed domains are currently missing.

The goal of this paper is to present a spatial evaluation of the SnowPappus blowing snow simulation framework (Baron et al., 2024), considering the uncertainty of precipitation estimates through different data sources that exhibit contrasted spatial patterns. The evaluation is based on Pléiades snow height maps (Deschamps-Berger et al., 2020) and Sentinel-2 and Landsat 8 snow melt-out dates (Gascoin et al., 2019) and covers three consecutive snow seasons. We discuss the relative influence of precipitation sources and blowing snow implementation on the simulated snow cover variability and how their interactions can affect the evaluation of a distributed snow model. By using two contrasting horizontal-resolution simulations, we also analyze how the unresolved spatial scales of blowing snow can affect simulation results at a 250 m spatial scale. Finally, we emphasize the main challenges to be solved for more advanced evaluations of spatialized snow simulations.

# 2 Data and methods

## 2.1 Study area

The simulation study site is located east of Grenoble in the French Alps (Fig. 1). It covers 902 km² including the Grandes Rousses and Arves massifs. This area exhibits a complex topography, with elevations ranging from 700 up to 3900 m a.s.l.; a wide range of snow and temperature conditions; and different types of landscape features such as valleys, forests, alpine pastures, lakes, and glaciers. In this area, the local knowledge of precipitation flow patterns tells us that most winter storms come from northwestern flows usually giving increased precipitation on western slopes. The study area includes the Col du Lac Blanc observatory (Guyomarc'h et al., 2019), where Baron et al. (2024) evaluated the blowing snow fluxes and occurrence simulated by SnowPappus. It counts 41 distinct ice patches or glaciers of various sizes as defined in the RGI Consortium (2017), forested areas, and three hydroelectric dams which underline the hydrological importance of the area.

As our simulation system is only able to represent snow in open areas, the forests, glaciers, lakes, and rivers in the study zone are masked in our simulations and observation datasets. For the forest mask, the BD FORET® V2 dataset has been used (IGN©, 2024a) with a masking threshold of 25 % of forested sub-pixel area. Waterways are masked following data from BD TOPO® (IGN©, 2024b). In the valley, the urban areas around the city of Bourg d'Oisans, Allemont,

and Saint-Michel-de-Maurienne are also discarded from the analysis.

## 2.2 Grid generation and spatial resolution

The simulations are based on two grid resolutions, 30 and 250 m. Both simulation grids are built using as reference the French 5 m RGE ALTI® digital elevation model (DEM) from IGN© (2021c). The 5 m high-resolution DEM is resampled using the average method to 30 and 250 m horizontal resolutions using GDAL/OGR contributors (2023). At 250 m resolution, the 902 km² full simulation area is composed of 14 443 grid points. At 30 m horizontal resolution, the test zone is composed of 1 005 699 simulation points.

A geomorphon classification (or landform) is performed on the 250 m resolution DEM for a more detailed analysis of simulation domain features. Geomorphons introduced by Jasiewicz and Stepinski (2013) are defined as the fundamental structural elements of a landscape. Here, we used the Whitebox Geospatial Inc. (2024) classification tool, an open-source software from Lindsay (2014). The algorithm is based on a line-of-sight analysis with similarities to the more classical topographic position index (TPI). This method identifies a set of topographic patterns corresponding to specific terrain attributes and landform types. The advantage of this type of classification is that it allows topographic information to be conveyed mimicking the result of a classification process carried out by a human analyst. Another advantage of the geomorphon classification is that it adapts to the surrounding terrain and can lead to the identification of landform elements regardless of their scale. Jasiewicz and Stepinski (2013) identify the 10 most common landform elements used for geomorphon classification (peak (summit), ridge, shoulder, spur (convex), slope, hollow (concave), footslope, valley, pit (depression), and flat; an illustration can be found in Fig. 3 of Jasiewicz and Stepinski, 2013). This allows us to generate a simple, intuitive, and scale-independent landform map of our simulation domain (Jasiewicz and Stepinski, 2013). Figure 2 illustrates the geomorphon classification result (called landform classification in the following) of each 250 m grid cell of our simulation domain in terms of the 10 most common landform elements. Figure 3 illustrates the elevation distribution of the landform elements of our simulation domain. No landform feature corresponding to Jasiewicz and Stepinski (2013) "shoulder", "footslope", and "flat" are found in our simulation domains. Finally, Appendix Fig. A1 demonstrates the good representativeness of the landform feature frequency of our domain compared with the entire French Alps and Pyrenees.

For our analysis, each 250 m simulation grid cell is classified in terms of landform types (see Fig. 2) and into 300 m elevation bands (see Fig. 1).
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Map of the simulation domain (red) and evaluation areas (Pléiades green and blue, snow melt-out date in red). Points of interest are shown as well as forest, glacier, and lake masks. The classification of 250 m pixels is visible with a 300 m step elevation in shades of blue. The summer aerial photography base map is from IGN© (2022).</figcaption></figure>

## 2.3 Crocus snow model

This paper uses the well-established 1D snow model Crocus (Brun et al., 1989; Vionnet et al., 2012), used operationally by the French weather forecast agency (Météo-France) in support of avalanche hazard forecasting (Morin et al., 2020). This snow model relies on a dynamical layering of the snowpack and an explicit representation of snow metamorphism and microstructure up to 50 snow layers (Vionnet et al., 2012). The Crocus snow model is coupled to the ISBA ground and land cover model and is embedded in the SURFEX modeling framework (Decharme et al., 2011; Masson et al., 2013). The snowpack is simulated on a squared grid of 250 m or 30 m spacing. Simulations are run over 3 years, starting from 2 August 2017 at 06:00 CET to 2 June 2020 at 06:00 CET, using a 15 min time step. The initial ground conditions are taken from a 10-year (2007–2017) simulation driven by the SAFRAN reanalysis (Vernay et al., 2022).

## 2.4 SnowPappus blowing snow model

The SnowPappus blowing snow model has been specially developed to describe snow variability due to blowing snow when the Crocus model is applied in the range of 25 to 250 m horizontal resolutions (Baron et al., 2024). SnowPappus model simulates blowing snow occurrence, horizontal

transport fluxes, and blowing snow sublimation rate on each grid cell as a function of 2D atmospheric forcing and snow surface properties. Then, a blowing snow mass balance is used to quantify eroded or accumulated snow amounts and modify the simulated snow profiles accordingly. The SnowPappus methods and parameterisations used to represent the different blowing snow processes are described in detail and compared to existing literature in Baron et al. (2024).

SnowPappus can run over multiple snow seasons and very large domains (about $10^5$ km$^2$ at 250 m resolution or $10^6$ simulation points), in a reasonable computing time (less than half a day on a single computing node, depending on simulation outputs).

## 2.5 Meteorological forcings

### 2.5.1 SAFRAN simulations

As most spatialized applications of the Crocus snow model over the French Alps, our simulations are forced by the SAFRAN meteorological reanalysis (Vernay et al., 2022). The SAFRAN analysis system uses as a background vertical profiles of air temperature, humidity, and the precipitation fields of the ARPEGE NWP system and assimilates near-surface meteorological observations to generate, with an hourly time step, the typical sets of inputs required by
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Map showing the result of the landform classification performed on our 250 m simulation domain. Each grid cell is classified in terms of the 10 most common landforms described by Jasiewicz and Stepinski (2013): peak (summit), ridge, shoulder, spur (convex), slope, hollow (concave), footslope, valley, pit (depression), and flat. In our domain, “flat”, “footslope”, and “shoulder” landforms are not present. The summer aerial photography base map is from IGN© (2022).</figcaption></figure>

a snow model such as Crocus (Durand et al., 1993). The SAFRAN reanalysis provides these meteorological inputs over 23 climatologically homogeneous areas (named massifs) covering the whole French Alps (Vernay et al., 2022). Here, only the SAFRAN reanalysis corresponding to the Grandes Rousses massif was used to build the simulation meteorological forcing. For a distributed use, the SAFRAN analysis is linearly interpolated as a function of elevation for each simulation point (as in Vionnet et al., 2016) to match the gridded geometry of our simulation domain. Solar radiations are projected according to the slope inclination and aspect of each pixel and masked in the case of shadows from the surrounding topography similarly to Revuelto et al. (2018) and Deschamps-Berger et al. (2022). In SAFRAN, all meteorological variables are assumed to be constant inside a massif for a given elevation. Whereas large-scale biases are low thanks to the assimilation of numerous gauge observations (Vernay et al., 2022), a large part of the spatial variability in precipitation remains unresolved under this assumption (Vionnet et al., 2016; Quéno et al., 2016; Vionnet et al., 2019; Deschamps-Berger et al., 2022). Therefore, other precipitation estimates are also considered in this study from the AROME NWP system and the ANTILOPE radar-based analysis. An illustration of the mean daily SAFRAN precip-

itation spatial variability can be seen in Appendix Fig. B1a. This map informs us of the spatial and inter-model variability and of the absolute precipitation values of each precipitation dataset.

## 2.5.2 AROME precipitation simulations

AROME is the high-resolution non-hydrostatic NWP system operated by Météo-France for short-range forecasts. It uses a 1.3 km grid spacing (Seity et al., 2011). Due to high temperature and radiative biases of AROME in mountain areas (Quéno et al., 2020; Gouttevin et al., 2023), and to focus our study on the uncertainty of precipitation spatial patterns rather than on all meteorological forcing uncertainties, we only consider the total precipitation amount (liquid and solid) from AROME in our simulations. Thus, the precipitation phase is set to solid if the SAFRAN air temperature is lower than 274.15 K and liquid if the air temperature is above 274.15 K. Although more advanced and continuous phase functions are available in the literature (Froidurot et al., 2014; Vionnet et al., 2022), this choice was made for consistency between all simulations, as this threshold is also used by the SAFRAN system. SAFRAN analysis data are used for all the other forcing variables. An illustration of the mean daily AROME precipitations can be seen in Appendix Fig. B1d.
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Distribution of simulation pixels in terms of 300 m step elevation bands and the landform classification. The number of pixels corresponding to each class is shown.</figcaption></figure>

### 2.5.3 ANTILOPE precipitation simulations

ANTILOPE is an hourly precipitation estimation analysis product (rain rate), combining radar and rain/snow gauges with 0.01° (1 km) resolution (Champeaux et al., 2009). ANTILOPE precipitation fields are the kriging result of available gauge measurements, using radar-estimated precipitation fields as external drift. In a mountainous context, the skill of the product can be affected by a partial or total mask of the radar beams and other common radar measurement errors (Faure, 2017; Yu et al., 2018), although the use of observations improves the precipitation estimates compared to raw radar precipitation. Since ANTILOPE does not discriminate solid from liquid precipitations, we used information from SAFRAN to classify rainfall and snowfall as well as all the other forcing variables. The precipitation phase is set to solid if the SAFRAN air temperature is lower than 274.15 K and liquid if the air temperature is above 274.15 K. SAFRAN analysis data are used for all the other forcing variables. An

illustration of the mean daily ANTILOPE precipitations can be seen in Appendix Fig. B1c.

## 2.6 DEVINE wind downscaling model

Wind fields are important drivers of turbulent and mass exchanges at the surface of the snowpack. Wind fields also profoundly determine the onset and evolution of drifting and blowing snow episodes. The interaction between mesoscale winds and local topography induces strong modifications in terms of both speed and direction that deeply influence the spatial variability in wind fields at the local mountain scale (Whiteman, 2000). Neither interpolated NWP systems operating with a kilometer resolution (Seity et al., 2011) nor massif-scale SAFRAN reanalyses (Vernay et al., 2022) would permit us to take into account the local influence of terrain on wind fields. As the local-scale wind patterns are crucial to determine high-resolution snow patterns (Musselman et al., 2015), we downscaled AROME wind fields (1.3 km grid spacing) to a 30 m grid using the DEVINE
downscaling model (Le Toumelin et al., 2022). This method leverages a convolutional neural network to emulate the behaviour of the complex atmospheric model ARPS (Advanced Regional Prediction System), previously run on a large number of synthetic topographies and in turn increase the spatial resolution of wind fields in complex terrain.

The DEVINE downscaling model has been previously evaluated over a large number of in situ stations (Le Toumelin et al., 2022, 2024), including stations located in our study domain. It has been shown that DEVINE wind outputs can represent wind acceleration along ridges and summits, deceleration windward, and some deflection around topographic obstacles, but it can not represent more complex wind patterns such as thermal winds or recirculation areas.

## 2.7 Model experiments

To disentangle the impacts of the simulation of blowing snow from the precipitation forcing, 250 m resolution simulations were run with and without the SnowPappus transport module and with the three different precipitation datasets, resulting in six different 250 m simulations. All simulations utilize the default options, including the GM98 parametrization of blowing snow occurrence, and the flux limiter is activated (Guyomarc'h et al., 1998; Baron et al., 2024). In addition, blowing snow sublimation is represented using parametrization adapted from the simplified blowing snow model (SBSM) (Essery et al., 1999; Baron et al., 2024). To investigate the role of unresolved transport variability in snow transport, two additional simulations were done at 30 m resolution with and without the SnowPappus transport module (see HR standing for high resolution). Due to the higher numerical cost of these simulations and the increased complexity of downscaling and validating weather forcing at this resolution, 30 m resolution simulations were only run with the SAFRAN precipitation fields.

As the 30 m resolution simulations are used in this paper to evaluate the contribution of unresolved spatial variability in 250 m simulations, we do not evaluate 30 m simulations at their native resolution but after a resampling on the 250 m resolution grid using the Earth System Modeling Framework (ESMF) first-order conservative method (Earth System Modeling Framework et al., 2023).

Regardless of the simulation configuration, 10-year SAFRAN spinup simulations (without SnowPappus) were performed to produce the ground's initial condition.

A summary of all the simulation experiments is found in Table 1.

## 2.8 Reference dataset and evaluation methods

### 2.8.1 Snow melt-out dates

We use the annual snow melt-out dates (SMODs) distributed by Theia (level 3 product of the snow collection (Gascoin

et al., 2019)). This product has a 20 m spatial resolution and indicates the last date of the longest continuous snow period. It is obtained by the linear interpolation in the time dimension of all the single-date (level 2) snow cover area products between 1 September and 31 August (Gascoin et al., 2019). These products are generated from Sentinel-2 and Landsat 8 images. The SMOD is given in days starting from 1 September. We resampled the SMOD to the 250 m resolution grid of our model using the Earth System Modeling Framework et al. (2023) first-order conservative method.

We use the same definition to compute the SMOD from the model output. To compute the longest continuous snow period, the snow height is first averaged by day and then a snow height threshold of 20 cm is applied, as this value was found to be optimal in the sensitivity analysis of Deschamps-Berger et al. (2022) at the same horizontal resolution. For this threshold, the accuracy is found to be 91.5% (See Fig. 10 of Gascoin et al., 2019). The forest, glaciers, and lakes on the domain are masked for both the observations and simulations using the mask defined in Fig. 1 and Sect. 2.1.

### 2.8.2 The Pléiades stereo imagery

In order to evaluate more directly the spatial variability in the simulated snow height, we used snow height maps derived from Pléiades stereo images (Deschamps-Berger et al., 2020; Marti et al., 2016). In both cases, a snow height is defined following a vertical axis. The Pléiades raw images have a 0.5 m horizontal resolution and are acquired on request. Snow height maps are obtained by combining two surface elevations of the same area, one with snow and one snow-free. Each surface elevation map is derived from pairs or triplets of Pléiades images. The reference snow-off surface elevation is calculated from Pléiades stereo images acquired on 28 September 2016 for the Galibier area and 24 October 2018 for the Grandes Rousses area. A winter image acquisition gives a second elevation map with snow on the ground. The difference between the two elevation fields provides a snow height map. The generated snow height map has a horizontal resolution of 2 m. Two observation zones have been selected based on available images. The first one covers 66 km² around the Col du Lac Blanc experimental site in the Grandes Rousses mountains. The second observation zone covers 168 km² around the Lautaret experimental site in the Galibier mountains. Two dates with snow on the ground have been analyzed for each observation zone. For the Grandes Rousses area, two images are available for a different year but at the same period of the season (Pléiades 1 on 13 May 2019 and Pléiades 2 on 4 May 2020). For the Galibier area, the two images give a view at different times of the same season (Pléiades 3 on 23 January 2018 and Pléiades 4 on 16 March 2018).

Due to the method used to reconstruct the snow height map, usual Pléiades observations contain a few no-data pixels where the quality of the observation is lower due to cliffs,
<table><caption>**Table 1.** Overview of the different simulation experiments, their model configurations, and naming. It is important to note that the computing resolution is different from the evaluation resolution. All simulations are evaluated at 250 m resolution (see Sect. 2.8).</caption><thead><tr><th>Experiment name</th><th>Precipitation forcing</th><th>Blowing snow transport mode</th><th>Computing resolution</th></tr></thead><tbody><tr><td>SAFRAN</td><td>SAFRAN</td><td>No transport</td><td>250 m</td></tr><tr><td>SAFRAN with transport</td><td>SAFRAN</td><td>Transport</td><td>250 m</td></tr><tr><td>SAFRAN HR</td><td>SAFRAN</td><td>No transport</td><td>30 m</td></tr><tr><td>SAFRAN HR with transport</td><td>SAFRAN</td><td>Transport</td><td>30 m</td></tr><tr><td>AROME</td><td>AROME</td><td>No transport</td><td>250 m</td></tr><tr><td>AROME with transport</td><td>AROME</td><td>Transport</td><td>250 m</td></tr><tr><td>ANTILOPE</td><td>ANTILOPE</td><td>No transport</td><td>250 m</td></tr><tr><td>ANTILOPE with transport</td><td>ANTILOPE</td><td>Transport</td><td>250 m</td></tr></tbody></table>

steep slopes, clouds, or shaded areas. The snow height map was filtered to exclude values out of the [-0.5 m, 20 m] range. According to previous Pléiades snow map evaluations (Deschamps-Berger et al., 2020), averaging 2 m Pléiades observation resolution to 250 m reduces the vertical snow height standard error to approximately 0.3 and 0.4 m. Negative snow heights are kept as removing them would slightly positively bias the averages. For comparison with snowpack simulations, the 2 m resolution Pléiades snow height map is resampled to the 250 m horizontal-resolution grid using a conservative method, similarly to Deschamps-Berger et al. (2020). The entire 250 m pixel is set to “no data” if more than 70 % of the 2 m pixels is “no data”; this threshold ensures the minimum representativeness of the pixel’s snow height.

The forest, glaciers, and lakes on the domain are masked for both the Pléiades observations and simulations using the mask defined in Sect. 2.1. An additional mask of the “no data” pixels in the Pléiades data is applied to the simulations for accurate comparison.

### 2.8.3 Observation data summary

All observations are summarized in Table 2. For the comparison with simulations, the observed SMOD and snow height are grouped by altitude in various elevation bands. The reference altitude used for the observation is the 250 m DEM described in Sect. 2.2. This elevation band classification is one of the most often used in the literature because the elevation is one of the main topographic drivers of snow cover variability (Grünewald et al., 2014; Vionnet et al., 2022; Monteiro et al., 2022; Tong et al., 2009; Deschamps-Berger et al., 2022; Vionnet et al., 2017). However, for a more advanced analysis of the snow cover spatial structure, we also present our results grouped according to the landforms presented in Sect. 2.2.

## 2.9 Synoptic scores

To summarize the simulation and observation comparison, a combination of three scores has been chosen. The mean

bias, standard deviation ratio, and spatial probability score (SPS) value are used to quantify the similarity between the simulated and observed snow height distribution.

The mean bias (MB) defined in Eq. (1) describes the mean bias between the simulation and the observation over a given area.

$$MB(S) = \langle S \rangle - \langle \emptyset \rangle, \quad (1)$$

where $\langle S \rangle$ represents the mean value of the simulation distribution and $\langle \emptyset \rangle$ the observation distribution.

The standard deviation ratio (or $\sigma$ ratio in the following) (Eq. 2) is the ratio between the simulated standard deviation and the observed one.

$$\sigma \text{ ratio}(S) = \frac{\sigma(S)}{\sigma(\emptyset)}, \quad (2)$$

with $\sigma(S)$ the simulated standard deviation and $\sigma(\emptyset)$ the observed standard deviation.

The spatial probability score (SPS) is a distance metric between two cumulative distribution functions (CDFs) defined in Eq. (3) as the quadratic discrepancy measure between the simulated CDF, denoted $F$, and $P_{\text{obs}}$ the empirical observed CDF.

$$SPS(F) = \int_{\mathcal{S}} (F(x) - P_{\text{obs}}(x))^2 dS, \quad (3)$$

with $x \in S$ and $S$ the volume of possible values of the indicator of interest. This score is derived from the continuous ranked probability score (CRPS) primarily used to evaluate probabilistic forecasts in ensemble forecasting (Candille and Talagrand, 2005). Here, it is used to compare two spatial distributions, similarly to the description in Goessling and Jung (2018), although using a non-binary probability field for the perfect observation. In short, this metric quantifies the similarity of two probability densities. It is close to the Wasserstein distance (Rüschendorf, 1985) used in Vionnet et al. (2021) although here of order 2. A low SPS value means a small distance between the observed and the simulated distributions. Two identical distributions give an SPS score of 0.
<table><caption>**Table 2.** Summary of satellite observations used in this study. The observed area, variables, resolutions, and dates are detailed for each observation. The aggregated resolution is a reminder of the resampling operation done on the observation data before comparison with simulations.</caption><thead><tr><th>Satellite observation</th><th>Observed area</th><th>Variable</th><th>Observation date</th><th>Native resolution</th><th>Aggregated resolution</th></tr></thead><tbody><tr><td>Pléiades 1</td><td>Grandes Rousses region (GR)</td><td>Snow height (m)</td><td>13 May 2019</td><td>2 m</td><td>250 m</td></tr><tr><td>Pléiades 2</td><td>Grandes Rousses region (GR)</td><td>Snow height (m)</td><td>4 May 2020</td><td>2 m</td><td>250 m</td></tr><tr><td>Pléiades 3</td><td>Galibier region (Glb)</td><td>Snow height (m)</td><td>23 January 2018</td><td>2 m</td><td>250 m</td></tr><tr><td>Pléiades 4</td><td>Galibier region (Glb)</td><td>Snow height (m)</td><td>16 March 2018</td><td>2 m</td><td>250 m</td></tr><tr><td>Sentinel-2 A</td><td>Entire simulation domain</td><td>SMOD (days)</td><td>1 September 2017 to 1 September 2018</td><td>20 m</td><td>250 m</td></tr><tr><td>Sentinel-2 B</td><td>Entire simulation domain</td><td>SMOD (days)</td><td>1 September 2018 to 1 September 2019</td><td>20 m</td><td>250 m</td></tr></tbody></table>

The SPS unit is identical to the unit of the target variable (m for snow height, days for SMOD).

# 3 Results

## 3.1 Spatial variability in the snow height

This section compares observed and simulated snow height from the eight experiments described in Table 1 to the four Pléiades observations (Table 2).

An example of simulated and observed snow height variability can be found in Fig. 4. In this figure, snow height is simulated on 13 May 2019 for the SAFRAN, AROME, and ANTILOPE precipitation forcings and observed for the Grandes Rousses area on the same date. For each precipitation forcing, the panel shows from one side the impact of the SnowPappus blowing snow model on the simulated snow height, where we can see areas of increased and decreased snow height, most visible along the ridges and peaks, and in each left sub-maps the same simulation without snow transport. It can be observed that the primary areas of model snow accumulation are located in regions of high observed snow height. However, the opposite is not true, many areas with high observed snow accumulation are not present in the simulations. Additional maps for 13 May 2019 comparing resolutions can be found in Appendix Fig. C1.

### 3.1.1 Grandes Rousses area

In this section, we compare the two Pléiades snow height images of the Grandes Rousses area to simulations. For a quantitative evaluation of these snow height maps, we summarized simulations and observations in terms of spatial distributions. For each Pléiades subdomain, simulated snow heights are discretized using 300 m elevation bands. Each snow height distribution describes the internal snow height spatial variability within each elevation band. The associated spatial distributions obtained in the Grandes Rousses region (green box in Fig. 1) are compared in Fig. 5a and b for 13 May 2019 (same data as Fig. 4) and 4 May 2020. The

skill of these distributions is summarized through synthetic scores presented in Fig. 6a.

Figure 5a shows a major impact of the precipitation input on the simulated snow height distribution, regardless of the use of snow transport in simulations. Precipitation inputs impact the simulated snow height distribution's shape, variance, and mean value. Simulation with AROME precipitations leads to the highest snow height on average for all elevation bands and the two evaluation dates. This result is especially notable for Fig. 5a at the elevation band of 2500 to 2800 m on the 13 May 2019 evaluation. The AROME precipitation gives a mean snow height of about 3 m representing a positive bias compared to the Pléiades observation means (1.87 m), consistent with Fig. 6a. The standard deviation of the AROME simulations overestimates Pléiades' standard deviation at low elevations on both dates (Fig. 6b, e). It is more realistic between 2500 and 2800 m and dependent on snow transport at higher elevations.

Simulations with the ANTILOPE precipitation, on the other hand, have a negative snow height bias compared to the observation. For the 2500 to 2800 m elevation band, the ANTILOPE simulations give a mean snow height of only 1.30 m. This negative bias tends to increase with altitude as illustrated in Fig. 6a. The standard deviation of ANTILOPE simulations is closer to observations than AROME simulations. Similar results are found for 4 May 2020 (Figs. 5b, 6d-f). For the SAFRAN and SAFRAN HR simulations, results are consistent between the simulation resolution (250 and 30 m) but with varying skill between both evaluation years. At all elevations, the SAFRAN and SAFRAN HR simulations show a lower bias compared to observations than simulations forced by AROME and ANTILOPE. However, the SAFRAN mean value of 4 May 2020 (around 2 m) underestimates the observed snow height. A similar result is found at both resolutions and with or without snow transport. At most elevations, the standard deviation is similar to ANTILOPE forcing and close to the Pléiades observations.

The representation of snow transport is found to only impact snow height bias for the upper elevation band (with a reduction in the mean snow height). The fact that this impact is beneficial or detrimental depends on the choice of the
<figure><img src="image_5.png" /></figure>

**Figure 4.** Maps illustrating raw simulation outputs and Pléiades observed snow height on 13 May 2019 at 10:00 CET. (a) SAFRAN 250 m simulation with snow transport and without snow transport in the left sub-map. (b) AROME 250 m simulation with snow transport and without snow transport in the left sub-map. (c) ANTILOPE 250 m simulation with snow transport and without snow transport in the left sub-map. (d) The Pléiades 1 observed snow height resampled to 250 m. The absence of snow in simulations is indicated in gray. For (d) the summer aerial photography base map is from IGN© (2022).

meteorological forcing. Below this elevation band, the im-
pact of snow transport on the bias is found to be shallow to
non-existent. Figure 5a also shows that for the three upper
elevation bands (above 2500 m), SnowPappus (non-hatched
violin plots) increases the simulated standard deviation on
average by 63.7% between snow height distribution with
and without snow transport (hatched violin plots). Figure 5b
(above 2500 m) shows a similar standard deviation increase
of 67.7% for the next year when adding snow transport in
the simulations. The impact of snow transport on the simu-
lated snow height is noticeable for all precipitation inputs.
Figure 6b and e show that the increase in standard deviation
(above 2500 m) introduced by the transport module is some-
times overestimated comparatively to the observation (σ ra-
tio > 1), depending on the date and precipitation forcing. Be-
low 2500 m, the snow height distributions with and without
snow transport, for a given precipitation input, are broadly
similar.

Finally, SPS scores (Fig. 6c, f) suggest that model perfor-
mances are mostly modified by SnowPappus at the highest
pixels of our domain and that simulation improvements or
degradations depend on combinations of model setups (trans-
port/no transport and precipitation forcings).

3.1.2 Galibier area

The second snow height evaluation is done in the Galibier
mountain area (blue box in Fig. 1). It uses two snow height
maps derived from Pléiades stereo images acquired in Jan-
uary and March 2018, in which snow cover extends to lower
elevations than in the Grandes Rousses snow height maps.
On both dates, a positive bias appears for the three lowest el-
evation bands, regardless of the forcing and decreasing with
elevation (Fig. 8a, d). For instance, the Pléiades observation
provides a mean snow height of 1.25 m at the 1900–2200 m
elevation band on 23 January 2018, while the simulation’s
mean values are all above 1.74 m of snow. At these eleva-
tions, the standard deviation is underestimated by all simula-
tions (0.71 m for the reference Pléiades observation, 0.16 m
for SAFRAN, and 0.32 m for AROME), without any signif-
icant impact of the snow transport module (see standard de-
viation ratio in Fig. 8b, e).

For the two upper elevation bands, the main features are
consistent with the evaluations for the Grandes Rousses area:
on both dates, the AROME simulations exhibit the high-
est mean snow height corresponding to a positive bias com-
pared to observations, while the ANTILOPE simulations ex-
hibit the lowest snow height corresponding to a negative bias.
Contrary to the previous section, SAFRAN simulations also
<figure><img src="image_6.png" /><figcaption>**Figure 5.** Spatial distribution of snow height (m) for 300 m elevation bands in Pléiades observations (gray) and all simulation configurations over the Grandes Rousses mountains on 13 May 2019 (Pléiades 1) (a) and 4 May 2020 (Pléiades 2) (b). The different colours represent the different precipitation forcing and horizontal resolutions, while hatching is used to distinguish between the activation and deactivation of the SnowPappus blowing snow scheme. Each violin contains a small inner black box representing the first and third quartile and a white dot representing the median distribution value. The number of samples for each subgroup is denoted by *N*.</figcaption></figure>

overestimate snow height in the Galibier area. For all precip-
itation forcing, the standard deviation of simulations with-
out transport (between 0.11 and 0.25 m in January at the up-
per level) is much smaller than in the observed spatial dis-
tribution (0.87 m). When the snow transport is activated, a
more realistic spatial variability is obtained, although it is
still underestimated below 2800 m elevation. This conclusion
remains valid on 16 March (Fig. 7b). As for the Grandes
Rousses evaluation, the simulations with the AROME pre-

cipitation have the highest spatial variability in snow height
at all elevations and for the two dates. Conversely, the sim-
ulations with SAFRAN have the smallest standard deviation
across elevation.

For the SPS score (Fig. 8c, f), the ANTILOPE forcing
produces very accurate scores (≤ 0.1 m) between 2200 and
2800 m and for both dates, in contrast to the results obtained
over the Grandes Rousses area. For the lowest elevation
band, AROME simulation with snow transport is the closest
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Scores comparing simulations with the Grandes Rousses Pléiades observations. Panels (a), (b), and (c) compare simulations with the 13 May 2019 Pléiades 1 image; panels (d), (e), and (f) are based on the 4 May 2020 Pléiades 2 image. (a), (d) Bias (m); (b), (e) standard deviation ratio (no unit); (c), (f) SPS (m). All scores are presented as a function of 300 m elevation intervals.</figcaption></figure>

to observations with an SPS score of 0.4 and 0.3 m in January and March. The simulations forced by SAFRAN precipitation provide the spatial distributions the furthest from Pléiades observations. This is an opposite conclusion to the results obtained for the Grandes Rousses region.

### 3.1.3 Summary: relative impacts of snow transport and precipitation forcing

When comparing the simulation performance between both simulation regions and all dates, common and different patterns arise. The impact of the simulation of snow transport on the bias is restricted to the upper altitude, close to moun-

tain summits. For all simulations, higher average snow height are obtained without transport than with transport. However, the added value of snow transport can not be established on these criteria as the bias highly depends on the precipitation forcing.

In all simulations, lower snow heights are obtained with ANTILOPE precipitation and higher snow heights with AROME precipitation, with increasing differences with elevation. However, upon examining the various bias scores (as shown in panels a and d of Figs. 6, 8, and 10), it is clear that no precipitation forcing produces a snow height simulation that is systematically less biased than the others over
<figure><img src="image_8.png" /><figcaption>**Figure 7.** Spatial distribution of snow height (m) for 300 m elevation bands in Pléiades observations (gray) and all simulation configurations over the Galibier on 23 January 2018 (Pléiades 3) (a) and 16 March 2018 (Pléiades 4) (b). Similar legend as in Fig. 5. Each violin contains a small inner black box representing the first and third quartile and a white dot representing the median distribution value. The number of samples for each subgroup is denoted by *N*.</figcaption></figure>

the different evaluation dates and areas. Across all elevation ranges, the mean driver of simulation bias appears to be the precipitation dataset choice.

The addition of snow transport increases the simulated standard deviation of snow height of 0.25 and reduces the

SPS value closer to the observation. The impact of transport on scores is mainly above 2500 m with a slope increasing with elevation. For these upper elevation bands, the snow transport prevails in the choice of the precipitation dataset on these criteria. For the lower elevation bands, the main driver
<figure><img src="image_9.png" /><figcaption>**Figure 8.** Scores comparing simulations with the Galibier Pléiades observations. Panels (a), (b), and (c) compare simulations with the 23 January 2018 Pléiades 3 image; panels (d), (e), and (f) are based on the 16 March 2018 Pléiades 4 image. (a), (d) Bias (m); (b), (e) standard deviation ratio (no unit); (c), (f) SPS (m). All scores are presented as a function of 300 m elevation intervals.</figcaption></figure>

of simulation standard deviation appears to be the choice of the precipitation dataset. The simulations using AROME precipitation have the biggest standard deviation value of all simulations.

### 3.1.4 Impact of subgrid variability

The effect of unresolved subgrid variability can be analyzed by comparing SAFRAN and SAFRAN HR simulations. Simulation maps illustrating the expected subgrid variability can be found in Appendix Fig. C1a, b, and c. At the 30 m grid resolution, the spatial variability in the topography (SAFRAN HR without transport vs. SAFRAN without transport) reduces the snow height and increases the 250 m scale spatial

variance at high elevations for the Grandes Rousses evaluations. The impact is smaller for the Galibier assessment. In Figs. 6 and 8, the difference in the different scores between a simulation with and without snow transport is smaller for simulations computed at a 30 m resolution, compared to the 250 m computed simulations (dashed and continuous orange lines are less distant or equivalent to the dashed and continuous green lines). Finally, it is observed that SAFRAN HR and SAFRAN simulations are relatively close in terms of snow height bias, $\sigma$ ratio, and SPS. Therefore, despite the significant impact of the simulation resolution at the process level, it is not obvious that the 250 m simulations would provide less accurate snow height spatial distributions than sim-
ulations computed at higher resolution and resampled on a
250 m grid.

## 3.2 Spatial variability in snow melt-out date (SMOD)

The evaluations of the SMOD allow us to extend the eval-
uations over our entire simulation domain. It also provides
more complementary information on the entire snow season
than the use of snow height snapshots. The associated spa-
tial SMOD distributions for the 2017–2018 and 2018–2019
years can be found in Fig. 9a and b. The synthetic summary
of simulation performances is given by Fig. 10.

3.2.1 Differences between years

There are some notable differences between the 2017–2018
and 2018–2019 snow season SMODs (Fig. 9). For the two
lower elevation bands, there are on average 15 d (Fig. 9a) and
20 d (Fig. 9b) SMOD differences between the 2 years for the
observation. A similar SMOD difference is found on average
for the simulations. However, for the 2018–2019 snow sea-
son (Fig. 9b), we see the observation and simulation distribu-
tions have a much bigger variance than the previous season.
Continuing the analysis, the observed median SMOD show
an earlier melt at all elevation for the 2018–2019 season. This
observation is even more pronounced for the three upper el-
vation bands. We also note that the distribution variability
is greater for this year. For the 2017–2018 snow season at
elevations between 2500 and 3400 m, there was an average
SMOD difference of 18 d between observations and simula-
tions. In the subsequent season, this difference increased to
an average of 39 d. In contrast, during the 2018–2019 snow
season, the simulations that included snow transport were,
on average, closer to the observed variance than for the other
seasons.

3.2.2 Impact of precipitation forcing

Above 2500 m, the simulation with transport and the ANTI-
LOPE precipitation has the largest distribution variance and
lower mean SMOD value. As a result, this is the only sim-
ulation that underestimates the mean SMOD in 2018–2019
(Fig. 10d) and overestimates its spatial variability (Fig. 10e).

For the other precipitation forcing, in 2017–2018 the SMOD is overestimated at low elevations, indicating that the snow melt-out occurred too late in the season compared to observation (Fig. 9b). Conversely, SMOD is underestimated at high elevations with a simulated melting up to 20 d before the observed one. In 2018–2019, SMOD was overestimated at all elevations. On average, the spatial variance of SMOD is underestimated except for simulations with the ANTILOPE forcing at the upper elevation where the σ ratio is found much above 1.

The SPS values of SMOD (Fig. 10c, f) finally provided
the opposite conclusion between both years, with the simu-
lations forced by AROME being the closest to the observed

spatial distribution at all elevations in 2017–2018 but the fur-
thest in 2018–2019 and conversely the simulations forced by
ANTILOPE the furthest to observations in 2017–2018 and
the closest in 2018–2019.

Analysis of wind speed and direction over the 3 simula-
tion years indicates that the values remained consistent. Ex-
amination of the yearly solid precipitation of the SAFRAN
forcing reveals that the 2018–2019 season had 21% less
snowfall compared to 2017–2018. This result is consistent
with a smaller amount of solid precipitation leading to a re-
duced amount of snow accumulation and the observed earlier
SMOD for the 2018–2019 snow season.

3.2.3 Impact of snow transport

Figures 9 and 10b and e show an increase in SMOD standard
deviation when activating snow transport in simulations, con-
sistent with the increase in snow height variance shown in
Sect. 3.1. The standard deviation of simulated SMOD with
transport is found to increase with altitude. This increas-
ing spatial variance of SMOD with altitude is also observed
with Sentinel-2 imagery, whereas simulation without trans-
port keeps a standard deviation rather similar. Nevertheless,
the impact of precipitation forcing on SMOD spatial vari-
ability is at least as high as the impact of snow transport.
Therefore, the added value of snow transport is still diffi-
cult to establish on those criteria, although simulations based
on SAFRAN and AROME forcing would suggest that the
SMOD spatial variability is improved by the snow transport
module.

The addition of snow transport also decreases the SMOD for all 250 m resolution simulations, with an increasing impact on altitude. Therefore, concerning the SMOD bias, the impact of snow transport can be either beneficial in 2018–2019 or detrimental in 2017–2018 depending on the sign of the general bias of the simulation without transport.

Except with ANTILOPE in 2017–2018, simulating snow
transport always results in improved SPS of SMOD. The
detrimental impact of snow transport at high elevations in
2017–2018 for ANTILOPE forcing is probably linked with
the much smaller estimates of precipitation by ANTILOPE
(Appendix Fig. B1), for which the reasons are discussed in
Sect. 4.2.

3.2.4 Impact of subgrid variability

The addition of snow transport in the 30 m resolution
SAFRAN HR simulations increases the SMOD for the 2
evaluation years, while in contrast snow transport decreases
the SMOD in the 250 m resolution SAFRAN simulation. In
terms of spatial variance, the increase in spatial variance by
snow transport is higher at 250 m resolution than at 30 m res-
olution. The direct impact of resolving a finer topography
in simulations without transport (differences between dashed
<figure><img src="image_10.png" /><figcaption>**Figure 9.** Snow melt-out date (SMOD) depending on elevation and simulation configuration. The simulations are compared with SMOD observations of the full simulation domain for the 2017–2018 (Sentinel-2 A) and 2018–2019 years (Sentinel-2 B). Each violin contains a small inner black box representing the first and third quartile and a white dot representing the median distribution value. The number of samples for each subgroup is denoted by *N*.</figcaption></figure>

orange and green lines) is less significant than on snow height evaluations.

## 3.3 Impact of landform types

In this section, the spatial simulated distributions of snow heights and SMOD are compared with satellite observations while grouped by landform types instead of elevation bands, as described in Sect. 2.2 and Fig. 2.

Figure 11 compares the January and March 2018 Galibier Pléiades images to snow height simulations grouped by landform types. Appendix Fig. D1 presents the same comparison using 13 May 2019 and 4 May 2020 Grandes Rousses Pléiades images. The comparison of simulation with and without transport shows that most of the snow height bias difference (Fig. 11a, d) is found on pixels classified as “peak (summit)” or “ridges”. However, this impact on bias can be either beneficial or detrimental depending on the bias of the simulation without transport, highly dependent on the precipitation forc-
<figure><img src="image_11.png" /><figcaption>**Figure 10.** Scores comparing simulated SMOD with Sentinel-2 imagery. (a), (b), (c) 2017–2018 year (Sentinel-2 A); (d), (e), (f) 2018–2019 year (Sentinel-2 B). (a), (d) Bias (days); (b), (e) standard deviation ratio (no unit); (c), (f) SPS (days). All scores are presented as a function of 300 m elevation intervals.</figcaption></figure>

ing, similar to the results described in the analyses based on elevation bands.

Simulations with snow transport show lower mean snow height in pixels classified as “peak (summit)” and “ridges” compared to simulations without transport. Conversely, simulations without snow transport exhibit lower snow height on “spur (convex)”, “slopes”, and “hollow (concave)”. This behaviour can be found in each 250 m simulation, regardless of the simulation forcing (Figs. 11, D1) except for the SAFRAN HR simulations.

It can be noticed that the ordering between the different precipitation inputs does not depend on the landform classes, while it was quite dependent on elevation bands.

For all precipitation datasets, the standard deviation of simulations without transport is much smaller than in the observed spatial distribution. The landforms where the activation of transport has the most impact on snow height variance are (by order of importance) “peak (summit)”, “ridges”, “spur (convex)”, “slopes”, and “hollow (concave)”. Compared to the observed variance, this behaviour is more realistic when snow transport is activated (Fig. 11b, e) for all precipitation forcings. While the simulated distribution variance better matches the order of magnitude of the observed one for “peak (summit)”, “ridges”, and “spur (convex)”, the simulated distribution variance of “slopes” and “hollow (concave)” is seen to be still half of the observed one. As for
<figure><img src="image_12.png" /><figcaption>**Figure 11.** Scores comparing simulations with the Galibier Pléiades observations. Panels (a), (b), and (c) compare simulations with the 23 January 2018 Pléiades 3 image; panels (d), (e), and (f) are based on the 16 March 2018 Pléiades 4 image. (a), (d) Bias (m); (b), (e) standard deviation ratio (no unit); (c), (f) SPS (m). All scores are presented as a function of landform pixel classifications.</figcaption></figure>

the evaluation using elevation bands, the simulations with the AROME precipitation have the highest spatial variability for the two dates and each landform type. Conversely, the simulations with SAFRAN and no snow transport have the smallest standard deviation across landforms.

The improvement in the SPS score (Fig. 11c, f) due to the representation of snow transport is major on “peaks (summits)” and “ridges” and noticeable on “spur (convex)” and “slope” classes. This improvement is obtained for all precipitation forcings (even in the ANTILOPE case where the bias has deteriorated) and is therefore mainly driven by the improvement in snow height variance. The SAFRAN HR simulations with and without transport have the lowest (better) SPS values for the two evaluation dates and most of the landform types, while it is not the least biased. Directly examining the SAFRAN HR cumulative distribution function re-

veals that the shape of the cumulative distribution is more similar to the observation despite being biased. This conclusion differs from the results obtained using a grid cell classification by elevation bands.

Figure 12 presents the 2-year evaluation using Sentinel-2 SMOD and landforms. Figure 12a and d show the impact of SnowPappus blowing snow simulations on simulated SMOD. As for snow height in Fig. 11, the impact on SMOD is most significant in areas classified as “peak (summit)” or “ridge”. The SMOD difference when using SnowPappus is found to go up to 20 d for “peak (summit)”. When comparing the 2 evaluation years, similar behaviour is found, except the 2018–2019 year appears uniformly more biased. Areas classified as “valley” and “pit” appear more biased for every precipitation forcing.
<figure><img src="image_13.png" /><figcaption>**Figure 12.** Scores comparing simulated SMOD with Sentinel-2 imagery. (a), (b), (c) 2017–2018 year (Sentinel-2 A); (d), (e), (f) 2018–2019 year (Sentinel-2 B). (a), (d) Bias (days); (b), (e) standard deviation ratio (no unit); (c), (f) SPS (days). All scores are presented as a function of landform pixel classification.</figcaption></figure>

In Fig. 12b and e, the impact of snow transport is again mainly present on “peak (summit)” and “ridges” landforms but with high variability between both years and a lower impact when AROME precipitation forcing is used (with already a higher spatial variability even without transport). In 2017–2018, the impact of snow transport is always beneficial for the SMOD σ ratio, as obtained for snow heights but for the second season, and it depends on the precipitation forcing. The simulated standard deviation of AROME simulations appears close to the observed one for all landforms except for pixels classified as “pit”. This better behaviour of AROME was not illustrated in the evaluation done using the elevation bands. The other precipitation forcing appears underdispersed with a standard deviation ratio of 0.5–0.75.

Looking at the SPS scores (Fig. 12c, f), the impact of snow transport is less clear than for snow heights, similar to the

results obtained by elevation bands. SPS values are smaller using landform classification than elevation bands. All simulations appear to have difficulties in representing SMOD distribution in pixels classified as “valley” and “pit” conversely to results of Fig. 11. The SPS scores are also more similar between the 2 evaluation years, contrary to the SPS score using the elevation band grouping. Figure 10c shows a distance increase between all simulations and observations distribution between 2500–3100 m. This effect is not found in Fig. 12c. In contrast to Fig. 11, the SAFRAN HR simulations have the overall highest SPS values, and AROME simulations’ σ ratios are found to be more realistic.

The impact of the group axis on results can be quantified by comparing identical data grouped differently. Figure 13 presents the mean SPS values of the eight simulations when grouped by elevation bands and landform types.
<figure><img src="image_14.png" /><figcaption>**Figure 13.** Plots of average SPS values for the different simulations when grouped with elevation bands or landform types. This graph illustrates the impact of the grouping axis (elevation or landform types) on snow height and snow melt-out date.</figcaption></figure>

We can see that SPS values are lower when grouped by landform types. We can see that grouping by elevation gives more homogeneous SPS scores for the different observation dates than grouping by landform where variation between dates is greater.

## 3.4 Impact of the SnowPappus model

Correlation between observations and the different simulations can be seen in Fig. 14. This figure illustrates scatter plots of snow height for the four different Pléiades observations for the eight different simulations (SAFRAN, SAFRAN HR, AROME, ANTILOPE, with and without activation of SnowPappus). For each simulation, we see more dispersed scatter plots for simulations with snow transport (greater variance). It can be noted that the mean difference in slope caused by the addition of snow transport is small (<0.02). Overall, AROME simulations have the steeper slopes. To quantify how correlated the simulation results are to the observations, we use the Pearson correlation coefficient. In Tables 3 and 4, we see Pearson correlation coefficients between Pléiades observations and simulations for each of the eight simulations. Correlation is computed two times for the entire area (all elevations) and restricted to elevation higher than 2700 m. Looking at the Pearson correlation for the entire Pléiades domain, the addition of snow transport in simulations decreases the correlation coefficient in 15 of 16 experiments. Restricting the correlation to pixels to elevation above 2700 m (pixels with high snow transport probability), the addition of snow transport in simulations increases the Pearson correlation for each experiment. This result is found statistically significant (with significance set at *p* value <0.05) in 15 out of 16 experiments. We note that for the entire domain, the variability in the correlation score for a single observation is greater according to the source of precipitation than for the addition of the snow transport process in the simula-

tions. To better disentangle the impact of precipitation forcing variability and the addition of snow transport, this result can be supplemented with Fig. 15, representing the mean SPS value of the six different 250 m simulations forced with SAFRAN, ANTILOPE, and AROME for elevations above 2700 m and the entire observation domain. In this figure, the addition of snow transport to the simulations clearly shows improvements in the mean SPS score on snow height and SMOD above 2700 m for all observations. SnowPappus lowers SPS values of 0.06 m and 0.85 d above 2700 m. The improvement in SPS score by adding snow transport is logically much smaller if we look at the total snow height and SMOD distributions without restricting the elevation. However, it is important to note that the addition of the SnowPappus model does not deteriorate the SPS scores.

# 4 Discussion

## 4.1 The added value of the SnowPappus model

In this study, we compare snow height maps simulated with and without transport using various sources for the precipitation forcing. The combined impact of precipitation forcing variability and wind-blown redistribution using the SnowPappus model is visible on snow height in Figs. 5 and 7 and SMOD in Fig. 9. For the higher altitudes, particularly around the summit pixels, the standard deviation of the simulated snow height and snow melt-out date distribution is primarily determined by blowing snow transport. The blowing snow transport impact is mostly restricted to “peak (summit)” and “ridges” landforms (Figs. 11, D1, 12) and more generally to the upper elevation bands (above 2700–2800 m). In addition, the snow height variance is seen to increase closer to the observation using SnowPappus’ simulated snow transport, regardless of the precipitation bias. Looking at the Pléiades bias values in panels a and d of Figs. 6 and 8, blowing snow lowers in value the mean bias of all simulations with snow transport compared to simulations without snow transport. This is probably a combined effect of increased snow density and mass loss from increased sublimation. This result is consistent with Vionnet et al. (2021). In Sect. 3.4, results show the SnowPappus model consistently improves the correlation between Pléiades observations and simulations above 2700 m. This improvement was found to be significant in 15 out of 16 experiments. In our simulation domain, the elevations of 2700 m and above correspond to areas where our snow transport model has a significant impact on the simulated snow height. For the Pléiades domain as a whole, the correlation coefficient is modestly reduced with the addition of snow transport. These results are a generalization of the correlation analysis done in Baron et al. (2024) and are consistent with the previous conclusions. However, the added value of SnowPappus can not be established only on these criteria as we find the simulated bias and standard deviation
<figure><img src="image_15.png" /><figcaption>**Figure 14.** Observed vs. simulated snow height scatter plots for the four different Pléiades observations (all elevations) and the different precipitation forcings. For each precipitation forcing, the simulation with snow transport is in blue and without transport in red. Linear regression line synthesizes distribution changes when adding snow transport. Additionally, the respective Pearson correlation coefficients ($\rho$) are displayed.</figcaption></figure>

dependent on the precipitation forcing. To better disentangle the effects of precipitation forcing variability and the addition of snow transport, the previous conclusions can be supplemented with Fig. 15, which shows the mean SPS value of the six different 250 m simulations forced with SAFRAN, ANTILOPE, and AROME. Based on these results, it can be concluded that the addition of the SnowPappus model to

our 250 m simulation improves the simulated snow distribution closer to the observed distribution for both snow height and SMOD. In the literature, Bernhardt et al. (2012) studied the influence of lateral snow redistribution using SnowTran-3D (Liston et al., 2007) and found an enhancement in the simulated snow spatial pattern using lateral snow redistribution. Using the Canadian Hydrological Model (Marsh et al.,
<table><caption>**Table 3.** Pearson correlation coefficients between observed and simulated snow height (Grandes Rousses). For each observation and simulation pair, a second correlation is computed only with pixel elevations (Z) above 2700 m. Underlined values represent *p* value >0.05. Bold values correspond to the best correlation value between simulation with and without the SnowPappus model.</caption><thead><tr><th>Forcing</th><th>Snow transport</th><th>Pléiades 1<br>13 May 2019</th><th>Pléiades 1<br>13 May 2019<br>Z>2700 m</th><th>Pléiades 2<br>4 May 2020</th><th>Pléiades 2<br>4 May 2020<br>Z>2700</th></tr></thead><tbody><tr><td>SAFRAN</td><td>Yes</td><td>0.744</td><td><strong>0.173</strong></td><td>0.766</td><td><strong>0.184</strong></td></tr><tr><td>SAFRAN</td><td>No</td><td><strong>0.751</strong></td><td><u>0.027</u></td><td><strong>0.795</strong></td><td><u>0.062</u></td></tr><tr><td>SAFRAN HR</td><td>Yes</td><td>0.739</td><td><strong>0.190</strong></td><td>0.803</td><td><strong>0.221</strong></td></tr><tr><td>SAFRAN HR</td><td>No</td><td><strong>0.741</strong></td><td><u>0.097</u></td><td><strong>0.820</strong></td><td>0.173</td></tr><tr><td>AROME</td><td>Yes</td><td>0.826</td><td><strong>0.169</strong></td><td>0.820</td><td><strong>0.211</strong></td></tr><tr><td>AROME</td><td>No</td><td><strong>0.835</strong></td><td><u>0.064</u></td><td><strong>0.850</strong></td><td>0.127</td></tr><tr><td>ANTILOPE</td><td>Yes</td><td><strong>0.592</strong></td><td><strong>0.174</strong></td><td>0.672</td><td><strong>0.198</strong></td></tr><tr><td>ANTILOPE</td><td>No</td><td>0.591</td><td><u>0.045</u></td><td><strong>0.693</strong></td><td>0.105</td></tr></tbody></table>

<table><caption>**Table 4.** Pearson correlation coefficients between observed and simulated snow height (Galibier). For each observation and simulation pair, a second correlation is computed only with pixel elevations (Z) above 2700 m. Underlined values represent *p* value >0.05. Bold values correspond to the best correlation value between simulation with and without the SnowPappus model.</caption><thead><tr><th>Forcing</th><th>Snow transport</th><th>Pléiades 3<br>23 January 2018</th><th>Pléiades 3<br>23 January 2018<br>Z>2700 m</th><th>Pléiades 4<br>16 March 2018</th><th>Pléiades 4<br>16 March 2018<br>Z>2700</th></tr></thead><tbody><tr><td>SAFRAN</td><td>Yes</td><td>0.492</td><td><strong>0.192</strong></td><td>0.586</td><td><strong>0.346</strong></td></tr><tr><td>SAFRAN</td><td>No</td><td><strong>0.599</strong></td><td><u>0.068</u></td><td><strong>0.675</strong></td><td>0.142</td></tr><tr><td>SAFRAN HR</td><td>Yes</td><td>0.535</td><td><u>0.110</u></td><td>0.626</td><td><strong>0.353</strong></td></tr><tr><td>SAFRAN HR</td><td>No</td><td><strong>0.598</strong></td><td><u>0.086</u></td><td><strong>0.664</strong></td><td><u>0.157</u></td></tr><tr><td>AROME</td><td>Yes</td><td>0.549</td><td><strong>0.189</strong></td><td>0.716</td><td><strong>0.414</strong></td></tr><tr><td>AROME</td><td>No</td><td><strong>0.600</strong></td><td><u>0.025</u></td><td><strong>0.744</strong></td><td>0.355</td></tr><tr><td>ANTILOPE</td><td>Yes</td><td>0.437</td><td><strong>0.181</strong></td><td>0.566</td><td><strong>0.422</strong></td></tr><tr><td>ANTILOPE</td><td>No</td><td><strong>0.509</strong></td><td><u>0.020</u></td><td><strong>0.606</strong></td><td>0.266</td></tr></tbody></table>

2020), Vionnet et al. (2021) find snowpack simulations without lateral snow redistribution unable to capture the spatial variability in snow cover in alpine terrain. In contrast, the addition of lateral snow redistribution is found to provide a better estimate of snow height across elevations and increase the simulated snow height distribution variance. Very recently, Quéno et al. (2023) evaluated the FSM2trans snow redistribution model (Quéno et al., 2023) and find the addition of snow redistribution necessary at hectometric or finer resolution to better represent snowpack heterogeneity. Snow heights simulated using lateral snow redistribution are found to have significantly better variability. Additionally, Quéno et al. (2023) conducted a model sensitivity analysis to horizontal resolution and established that, although the most realistic simulation pattern is found at 25 m resolution, the addition of lateral snow redistribution at coarser resolution is also having a positive impact on snow distribution and variance. Hence, we estimate the results obtained in this study

on the impact of the SnowPappus snow transport model are compatible with the literature.

## 4.2 Sensitivity to precipitation forcings

In our simulations, the mean and median snow height and SMOD values are primarily driven by the spatial variability in precipitation forcing, as shown in Figs. 5, 7, and 9. Precipitation forcing is also found to be the main driver of simulated snow height distribution standard deviation across low to medium elevations. As the Crocus snow model is known to produce realistic snow heights and snow water equivalent when forced by well-controlled meteorological forcing (Lafaysse et al., 2017; Menard et al., 2021), the large biases of snow height and SMOD observed in the snow simulations presented most probably reflect intrinsic biases of the different precipitation datasets. In panels a and d of Figs. 6 and 8, for a simulation with a precipitation dataset leading to an already negatively biased snow height (like ANTILOPE),
<figure><img src="image_16.png" /><figcaption>**Figure 15.** Plots of average SPS values for the three different 250 m forcings (SAFRAN, ANTILOPE, AROME) compared to observations. SPS values are computed for pixels with elevations above 2700 m (square) and pixels of the entire observation domain (all elevations, circle). This graph summarizes the impact of the SnowPappus transport model on snow height and snow melt-out date distributions.</figcaption></figure>

the addition of snow transport further degrades the overall snow height simulation bias. Similar behaviour is found for the SMOD in Fig. 10a and d with an exception for the 30 m simulation where it increases the absolute SMOD bias. Thus the beneficial or detrimental impact of blowing snow transport on mean snow height or SMOD bias is mainly dependent on the precipitation bias. This result is consistent with Raleigh et al. (2015), Schlögl et al. (2016), and Günther et al. (2019), also finding the prevalence of meteorological forcing errors in snow simulations. More specifically, when looking at the simulation forcing with the ANTILOPE precipitation source, we found that across all elevation bands, the ANTILOPE forcings have the lowest snow height and SMOD values. A negative bias trend is seen with elevation (Fig. B1). The ANTILOPE precipitation product is primarily based on radar reflectivity measurements, a technique known to be more challenging in mountainous terrains (Yu et al., 2018; Faure et al., 2019; Foresti et al., 2018; Germann et al., 2022). Meteorological radars are based on the backscatter measurement of several conic-shaped microwave electromagnetic impulses produced by the radar itself at several elevation angles (conic-shaped radar beams). When a mountain intercepts a radar beam, the reflected signal (called ground clutter) is unusable and rejected by the ANTILOPE algorithm. If the mountain intercepts less than 70% of the beam section, the signal reflected further away is simply corrected to account for the loss of information but is still used for precipitation rate estimation. Precipitation estimation over pixels affected by ground clutter thus comes from beams that can be several thousands of meters higher than the mountain, often above the base of precipitation clouds.

This leads to an underestimation of the precipitation rate over mountain ridges (Appendix Fig. B1). One possible compromise would consist of using the observed spatial variability from ANTILOPE while correcting the elevation-dependent precipitation bias.

In contrast to ANTILOPE, our simulations with the AROME precipitation source lead, on average, to the highest snow height across the elevation bands. It is particularly visible in Fig. 4b, where squared patterns can be seen and come from the 1.3 km AROME pixels forming distinctive precipitation patterns (see Appendix Fig. B1d). This behaviour of excessive solid precipitation at high elevations is known and coherent with other snow cover studies using AROME as the precipitation source (Quéno et al., 2016; Vionnet et al., 2019; Monteiro et al., 2022). Multifactorial causes were discussed by Monteiro et al. (2022) and Gouttevin et al. (2023). The variable biases between the western (Grandes Rousses) and eastern (Galibier) evaluation areas suggest that orographic enhancement of precipitation might be especially overestimated by AROME for the areas the most exposed to the prevailing westerly flows.

Finally, the two SAFRAN simulations show a good agreement in terms of SMOD and snow height biases but are spatially under-dispersed by design. The SAFRAN meteorological reanalysis (Vernay et al., 2022) uses one vertical profile per climatologically homogeneous area, interpolated to all simulation points with altitude (see Sect. 2.5.1). Within one of these areas, the spatial variability is thus limited to a vertical gradient as described in Vionnet et al. (2016). This also explains why opposite snow height biases are obtained with
this forcing between western (Grandes Rousses) and eastern (Galibier) areas.

## 4.3 Sensitivity to horizontal resolution

Numerical simulations rely on numerical discretization and parameterization. The choice of a spatial resolution inevitably introduces errors. This choice depends on a trade-off between scientific arguments, computing resources, and external factors. As described in Baron et al. (2024), the 250 m resolution choice of the SnowPappus scheme was motivated by the desire for applicability in large-extent systems (50 000 to 100 000 km², multi-decadal simulations). The numerical resolution choice impacts simulated slope inclinations and aspects and can lead to inaccuracy in the radiative balance (Baba et al., 2019). The snow radiative balance is particularly heterogeneous in mountainous areas due to complex topography and can lead to significant changes in snow properties and melting rates. Comparing simulations computed at 250 and 30 m resolution provides an analysis of the unresolved subgrid variability. A direct comparison between Appendix Fig. C1a, b, and c reveals the expected unresolved subgrid variability. Snow height patterns in the 250 m simulations (Appendix Fig. C1a) exhibit steeper variations than the 250 m re-gridded maps obtained from 30 m simulations (SAFRAN HR, Appendix Fig. C1b) where the snow patterns are found to be smoother. To understand the origin of this subgrid variability, Appendix Fig. C1c illustrates the raw 30 m simulation where snow patterns are more sophisticated. As expected (Baba et al., 2019), resolving the 30 m topography in simulations without snow transport increases the spatial variance at high elevations and reduces the 250 m averaged snow height (Sect. 3.1.4). Increasing the resolution in simulations with snow transport increases the SMOD (Sect. 3.2.4) but reduces the expected snow height variability gain compared to 250 m simulations. The complex behaviour of resolution can be summarized by comparing the mean SPS value across all elevations of the SAFRAN and SAFRAN HR simulations, as in Appendix Fig. E1. This figure shows a very small mean SPS change with resolution. The mean difference is found to be lower than 0.01 m for the Pléiades observations and lower than a day (0.98) for SMOD. Until now, the resolution-dependence of snow transport models was only documented at higher resolutions (Grünewald et al., 2013; Marsh et al., 2020; Mott et al., 2010; Schneiderbauer and Prokop, 2011), but results recently submitted suggest that the resolution-sensitivity of SnowPappus in this range of spatial scales is similar to what is obtained from another snow transport model at similar resolutions (Quéno et al., 2023). Quéno et al. (2016), Bellaire et al. (2014), and Schirmer and Jamieson (2015) also emphasize the relationship between the available resolution of the atmospheric forcing and the resulting snow simulation leading to the recommendation of a forcing spatial resolution of 1 km or less, although evaluations of these data are very challenging at this

scale. A much lower resolution is needed to capture the relationship between topography, wind field and snow accumulation (Dadic et al., 2010) or for a complete representation of wind-driven processes (Mott et al., 2008, 2014; Bernhardt et al., 2010; Vionnet et al., 2017).

## 4.4 Impact of the grid cells classification in spatial evaluations

Classically in the snow modeling community, the quantitative spatial snow evaluations are performed by grouping simulation points (Grünewald et al., 2014; Vionnet et al., 2022; Monteiro et al., 2022; Tong et al., 2009; Deschamps-Berger et al., 2022; Vionnet et al., 2017) because (1) it allows summarizing large simulation datasets with more concise diagnostics and (2) these evaluations are less demanding for models on large areas than pixel-to-pixel evaluations. In Sect. 3.3, results are presented using landform classification instead of elevation bands. The results are partly correlated because summits and ridges mainly cover the highest elevations of the domain (Fig. 3). Thus, we consistently obtain snow transport mainly affecting snow height and SMOD on the upper elevation bands and for the geomorphons “peaks (summit)” and “ridge”. However, some conclusions significantly differ between both space classifications, especially for the σ ratio score and the SPS. A notable difference is found in the SPS snow height comparison. SPS scores are found to show lower variability for the different precipitation forcings and overall smaller SPS values using landform groupings than elevation bands (Fig. 11, Appendix Fig. D1c, f). This result can be supplemented with Fig. 13 showing the mean SPS value being consistently better using landform groups. This suggests simulated and observed distributions are closer when grouped using landforms than elevation bands. Indeed, for the landform grouping, each different landform distribution includes a relatively large spatial variability (intra-class variance) mostly due to elevation, slope, aspect, and precipitation variability (Clark et al., 2011; Freudiger et al., 2017). On the other hand, using the elevation grouping the intra-class variance is mostly caused by other processes than elevation. This might suggest that our simulations better capture the altitudinal variability in snow than the other processes and the topographic-dependent variability. This result clearly illustrates the sensitivity to the choice of grid cell classification groups of the scores used. The landform classification also provides complementary insights into the behaviour of the SnowPappus model. In Fig. 11 and Appendix Fig. D1, SnowPappus reduces the mean snow height on pixels classified as “peak (summit)” and “ridges” and increases the snow height on “spur (convex)”, “slope”, and “hollow (concave)” comparatively to simulations without snow transport. This behaviour is the result of the wind spatial patterns obtained from Le Toumelin et al. (2024). It can be understood as snow on “peak (summit)” and “ridges” (ablation areas) being transported to pixels classified as “spur (convex)”,
“slope” or “hollow (concave)” (deposition areas) due to respective wind acceleration and deceleration on these pixels. This behaviour is similar between all 250 m resolution simulations, while different patterns are obtained from the 30 m SAFRAN HR simulation. Indeed, the finer spatial scales resolved at 30 m resolution produce more complex spatial patterns of wind fields and snow transport. When averaged at 250 m resolution, these simulations produce a more subtle influence of snow transport than simulations at the coarser 250 m resolution, with a less direct influence of the 250 m topography on the simulated snow patterns.

## 4.5 Limitations of simulations

In this study, we tested the sensitivity of simulated spatial snow height and SMOD variability to the variability in precipitation, snow transport, and spatial resolution, but numerous other processes shape snow spatial variability in the field, as reviewed in Mott et al. (2018). Those processes occur at different scales, going from the regional scale to the slope and lower scale. It is worth noticing that neither preferential deposition nor the seeder-feeder mechanism is represented in our simulations. Additional errors in snow-wind interaction come from process modeling. In our simulation, only the blowing snow redistribution process is represented, but the preferential deposition and snowfall enhancement are not. Blowing snow transport parameterization also suffers from uncertainties. The parameterizations used in our model were developed for flat areas and the question remains on using them in complex turbulent mountain topography. As discussed by Aksamit and Pomeroy (2018), turbulence, wind gusts, and eddies are some of the main small-scale contributors to blowing snow. In any case, the DEVINE downscaled wind does not yet take into account wind turbulence and recirculation zone. Further studies are needed to quantify the relative impact of recirculation zones over lee deceleration at the 250 m spatial resolution used in this work. Moreover, in our simulations not all cover types and processes are represented. Forests and glaciers are masked out for comparisons in addition to unusable observed pixels (due to topographic shadows or other observation limitations). The necessity of using these masks inevitably has an impact on the results by reducing the number of pixels in use for the evaluation. For masked forest pixels, we know that snow transport is strongly inhibited in forested areas leading to the creation of snow accumulation at the forest edges (Bernhardt et al., 2012). It is not possible to account for this phenomenon solely through post-processing masks. In our simulation domain, the main forested areas are located at low elevations where snow transport is minimal or non-existent. On the other hand, glaciers are generally located in areas prone to blowing snow. Removing glaciers from the analysis significantly decreases the total number of pixels subject to snow transport at higher altitudes, but it does not introduce bias to pixels at the edges of the mask. Although the simulated and observed snow heights

are not directly usable, it is important to note that snow transport is still simulated on masked glacier pixels and thus can contribute to neighbour pixels' snow balance.

## 4.6 Uncertainties and limitations of observations

In this work, simulations are evaluated against spatialized snow observations, giving a larger perspective than point-scale snow observations. Indeed, although point-scale snow observations are commonly used to evaluate spatialized simulations (Horton and Haegeli, 2022; Vernay et al., 2022; Vionnet et al., 2022; Mott et al., 2023; Marsh et al., 2020), their spatial distribution and representativeness is quite low in mountains (Grünewald et al., 2015; Pomeroy et al., 2009; Pepin et al., 2015), while satellite snow observations like Pléiades snow height and Sentinel-2 SMOD have a much better spatial coverage and representativeness. Nonetheless, some errors and limitations also come from the observation side. The Pléiades observations are on-demand snapshots of snow height over a predefined extent and at a specific date. By design, this limits the usage of Pléiades images and the extent to which our conclusions stay valid. Pléiades observations of snow height can be affected by shadows, steep slopes, and the reference image used in the processing method. These limitations restrict the domain of observations and make it difficult to retrieve accurate snow height on glaciers. The Pléiades horizontal resolution of 2 m causes challenges in comparison to coarser-resolution simulations. Unsimulated subgrid processes are captured by the observation, which complexifies the analysis (Fig. 4). For the Sentinel-2 observations, the major limitation comes from the snow detection algorithm. The presence or absence of snow is not directly comparable to our simulation output. In Sect. 2.8.1, the presence or absence of snow is expressed as a function based on a given simulated snow height threshold. This raises the question of the sensitivity to this threshold (Gascoin et al., 2019; Hofmeister et al., 2022).

## 5 Conclusions

We conducted a spatialized evaluation of the blowing snow model SnowPappus (Baron et al., 2024) joined with a sensitivity analysis to alpine precipitation variability on a 902 km² area representative of the French Alps and the Pyrenees in terms of landforms. We evaluated the simulated snow height and snow melt-out date (SMOD) using Pléiades and Sentinel-2 satellite products. The Pléiades comparisons were conducted in two different areas, on four different dates, and during three different snow seasons. Eight different snow simulations have been run over three snow seasons with three different precipitation forcings and two simulation resolutions (Table 1). We performed simulation analysis with the aim of disentangling simulation error contributions from the SnowPappus blowing snow model, the precipitation forcing
variability, and the unresolved subgrid variability. The main conclusions of this study are as follows:

- Simulations without snow transport are found to be unable to capture the spatial variability in snow cover in alpine terrain. The addition of SnowPappus' snow transport in the simulations results in a more physically realistic assessment of snow height and SMOD, regardless of precipitation forcing, and leads to simulated snow height and SMOD closer to the observations. Furthermore, the use of the SnowPappus model increases the variance of simulated snow height regardless of precipitation forcing bias.

- Precipitation errors are identified as the main source of bias and standard deviation in snow height and SMOD at low to medium altitudes. This result underlines that the greatest care must be taken to obtain the most accurate precipitation fields for snow simulation. It also suggests that assessments of spatial snow simulations that do not account for precipitation uncertainty are unlikely to provide informative insights regarding the accuracy of particular simulated snow processes.

- Blowing snow transport impact on snow height, SMOD bias, and standard deviation prevail for high altitudes or pixels classified by the landform classifier as “peak (summit)” and “ridges”.

- The changes in simulation snow height and SMOD mean SPS when increasing computation resolutions from 250 to 30 m are minimal. For Pléiades observations, the change is less than 0.01 m, and for SMOD, it is less than 1 d. These changes in spatial variance and bias are also lower than the simulated spatial variability due to precipitation inputs and the impact of blowing snow in the affected areas.

- The addition of the SnowPappus blowing snow model to a 250 m Crocus snowpack simulation consistently increases and improves the simulated snow spatial variability at high elevations in accordance with observations. However, improvements in snow height and SMOD biases are only obtained using precipitation forcing that does not suffer from strong negative biases.

The findings of this study show promising results of using the SnowPappus blowing snow model for large-scale modeling of alpine snowpack. Future work will focus on adding to the simulation workflow assimilation methods to improve precipitation estimates. Our results suggest that snow transport impact on snow height at 250 m resolution is significant at high elevations and therefore needs to be accounted for. However, unresolved processes at this spatial scale, parameterization uncertainty, and evaluation challenges are responsible for large uncertainties that will need to be appropriately

quantified in future works, especially with the goal of designing a snow ensemble assimilation framework (e.g., Cluzet et al., 2021; Deschamps-Berger et al., 2022).
Appendix A: Additional domain characteristics

<figure><img src="image_17.png" /><figcaption>**Figure A1.** Landform classification frequency (at 250 m) for our simulation domain and the entire French Alps and Pyrenees (as an area defined in Vernay et al., 2022). We found a good agreement between our simulation pixel classification and the frequency of the French Alps and Pyrenees. We note that using 250 m resolution, no pixel corresponding to shoulder, footslope, and flat is classified as such.</figcaption></figure>
Appendix B: Additional map of mean daily
precipitation variability

a)

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

b)

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

c)

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

d)

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

**Figure B1.** Map of the mean daily precipitation over 2017–2020 (in mm d⁻¹) for the four different forcings. This illustrates the spatial variability in each precipitation dataset: (a) SAFRAN, (b) SAFRAN HR, (c) ANTILOPE, (d) AROME. In each map, the mean solid fraction of precipitation is displayed as a contour plot (i.e., 1 = only solid precipitation). The summer aerial photography base map is from IGN© (2022).
Appendix C: Additional map of simulated snow height
for Pléiades 1, 13 May 2019

<figure><img src="image_22.png" /><figcaption>**Figure C1.** Map of simulated snow height on the 13 May 2019 at 10:00 CET. (a) SAFRAN 250 m simulation with snow transport and without snow transport in the left sub-map. (b) SAFRAN HR simulation (30 m computed resampled to 250 m) with snow transport and without snow transport in the left sub-map. (c) Raw SAFRAN 30 m simulation with snow transport and without snow transport in the left sub-map. The absence of snow is indicated in gray.</figcaption></figure>
Appendix D: Additional landform analysis

<figure><img src="image_23.png" /><figcaption>**Figure D1.** Scores comparing simulations with the Grandes Rousses Pléiades observations. Panels (a), (b), and (c) compare simulations with the 13 May 2019 Pléiades 1 image; (d), (e), and (f) are based on the 4 May 2020 Pléiades 2 image. (a), (d) Bias (m); (b), (e) σ ratio (no unit); (c)–(f) SPS (m). The important standard deviation ratio values found for the “pit” landform type in (b) and (e) are explained because almost all pixels of this type have no observed snow, whereas the simulations give snow on every pixel. This is particularly visible in panels (d) and (e) where ANTILOPE, SAFRAN, and SAFRAN HR have almost no simulated snow and a fair standard deviation ratio except for ANTILOPE. All scores are presented as a function of landform pixel classification.</figcaption></figure>
## Appendix E: Additional resolution analysis

<figure><img src="image_24.png" /><figcaption>**Figure E1.** Plots of average SPS values for the SAFRAN and SAFRAN HR simulations when grouped with 300 m elevation bands. This graph illustrates the impact of the computing horizontal resolution on snow height and snow melt-out date.</figcaption></figure>

**Code and data availability.** The SnowPappus blowing snow model is developed in the framework of the open-source SURFEX project. The source files of the SURFEX system (Crocus snow model, ISBA ground model, and SnowPappus model) are provided at https://doi.org/10.5281/zenodo.7687821 (Baron et al., 2023a) to guarantee the permanent reproducibility of results. However, we recommend that potential future users and developers access the code from its Git repository to benefit from all tools of code management (history management, bug fixes, documentation, interface for technical support, etc.). This requires a free and quick registration. The procedure is described at https://github.com/UMR-CNRM/snowtools/blob/master/doc/source/misc/surfex-install.rst (last access: 24 June 2024). The version used in this work is tagged as SnowPappus-v1.0. A user manual describing the SURFEX namelist options related to SnowPappus is available at https://doi.org/10.5281/zenodo.7687830 (Baron et al., 2023b). More general information about SURFEX use can be found at https://github.com/UMR-CNRM/snowtools (last access: 24 June 2024).

The code and model weights used for DEVINE are available for free at https://github.com/louisletoumelin/neural_network_and_devine/ (last access: 24 June 2024; DOI: https://doi.org/10.5281/zenodo.10594274, Le Toumelin, 2024).

The WhiteboxTools Open Core software used to compute the landform classification is available for free at https://www.whiteboxgeo.com/manual/wbt_book/available_tools/geomorphometric_analysis.html#Geomorphons (Whitebox Geospatial Inc., 2024).

The Earth System Modeling Framework (ESMF) software is available for free at https://earthsystemmodeling.org/ (last access: 24 June 2024; DOI: https://doi.org/10.5281/zenodo.11205527, Theurich et al., 2024).

The GDAL software is also available for free at https://doi.org/10.5281/zenodo.8340595 (Rouault et al., 2023).

The DEM used in this study originates from the IGN® website. Tiles can be downloaded for free at https://geoservices.ign.fr/rgealti (IGN©, 2021c).

Forest mask originates from IGN® and can be downloaded for free at https://geoservices.ign.fr/bdforet (IGN©, 2024a).

Waterways, lakes, and cities' spatial extent originate from IGN® and can be downloaded for free at https://geoservices.ign.fr/bdtopo (IGN©, 2024b).

The glacier masks can be downloaded for free from the Randolph Glacier Inventory at http://glims.colorado.edu/glacierdata/ (last access: 24 June 2024; DOI: https://doi.org/10.7265/N5V98602, GLIMS and NSIDC, 2005).

The SAFRAN reanalysis dataset is available for free on the AERIS data center via the following link https://doi.org/10.25326/37#v2020.2 (CNRM/Centre d'Etudes de la Neige, 2019).

The AROME and ANTILOPE precipitation dataset can be requested online https://donneespubliques.meteofrance.fr/ (last access: 24 June 2024).

The Theia level 3B snow products can be accessed at https://doi.org/10.24400/329360/F7Q52MNK (Theia Snow, 2019).

The raw Pléiades images are not publicly available. The 250 m averaged Pléiades images can be found at https://doi.org/10.5281/zenodo.10037253 (Haddjeri et al., 2023).

**Author contributions.** AH ran the simulations, collected and processed the different data, wrote the paper, and made the spatial evaluations and the scientific choices. AH and MB developed jointly the SnowPappus code with equal contribution. MB helped write the Crocus and SnowPappus sections. ML supervised the work, helped a lot with the scientific choices, and provided technical support on SURFEX/Crocus. ML and MD were extensively involved in the proofreading process and the scientific choices. LLT provided the DEVINE downscaled wind forcing and helped in writing the wind DEVINE section. MV provided the ANTILOPE precipitation forcing and helped the writing of the ANTILOPE section. CDB and SG provided the Pléiades images and the stereo imagery treatment. CDB, SG, and VV helped with the evaluation methodology, participated in theoretical discussions, and proofread the article.

**Competing interests.** At least one of the (co-)authors is a member of the editorial board of *The Cryosphere*. The peer-review process was guided by an independent editor, and the authors also have no other competing interests to declare.

**Disclaimer.** Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors.

**Acknowledgements.** This investigation was done using Météo-France's CNRM and CEN computational and human resources. CNRM and CEN are part of LabEx OSUG@2020. We gratefully acknowledge Sabine Radanovics, Diego Monteiro, and
Matthieu Fructus for their expertise. We thank the editor Masashi Niwano and the different reviewers for their careful review process, which significantly improved this paper.

*Financial support.* This research has been supported by the Region Auvergne Rhône-Alpes (France) SENSASS project, CNES TOSCA, through the TRISHNA-Cryosphere project.

*Review statement.* This paper was edited by Masashi Niwano and reviewed by Hiroyuki Hirashima and one anonymous referee.

## References

Aksamit, N. O. and Pomeroy, J. W.: Scale Interactions in Turbulence for Mountain Blowing Snow, J. Hydrometeorol., 19, 305–320, https://doi.org/10.1175/JHM-D-17-0179.1, 2018.

Amory, C., Kittel, C., Le Toumelin, L., Agosta, C., Delhasse, A., Favier, V., and Fettweis, X.: Performance of MAR (v3.11) in simulating the drifting-snow climate and surface mass balance of Adélie Land, East Antarctica, Geosci. Model Dev., 14, 3487–3510, https://doi.org/10.5194/gmd-14-3487-2021, 2021.

Anderson, B. T., McNamara, J. P., Marshall, H.-P., and Flores, A. N.: Insights into the physical processes controlling correlations between snow distribution and terrain properties, Water Resour. Res., 50, 4545–4563, https://doi.org/10.1002/2013WR013714, 2014.

Baba, M. W., Gascoin, S., Kinnard, C., Marchane, A., and Hanich, L.: Effect of Digital Elevation Model Resolution on the Simulation of the Snow Cover Evolution in the High Atlas, Water Resour. Res., 55, 5360–5378, https://doi.org/10.1029/2018WR023789, 2019.

Bales, R. C., Molotch, N. P., Painter, T. H., Dettinger, M. D., Rice, R., and Dozier, J.: Mountain hydrology of the western United States, Water Resour. Res., 42, W08432, https://doi.org/10.1029/2005WR004387, 2006.

Baron, M., Haddjeri, A., Lafaysse, M., Le Toumelin, L., Vionnet, V., and Fructus, M.: Supplementary to “SnowPappus v1.0, a blowing-snow model for large-scale applications of Crocus snow scheme”: SURFEX codes and dependancies, Zenodo [code], https://doi.org/10.5281/zenodo.7687821, 2023a.

Baron, M., Haddjeri, A., Lafaysse, M., Le Toumelin, L., Vionnet, V., and Fructus, M.: Supplementary material to “SnowPappus v1.0, a blowing-snow model for large-scale applications of Crocus snow scheme”, user manual, Zenodo, https://doi.org/10.5281/zenodo.7681340, 2023b.

Baron, M., Haddjeri, A., Lafaysse, M., Le Toumelin, L., Vionnet, V., and Fructus, M.: SnowPappus v1.0, a blowing-snow model for large-scale applications of the Crocus snow scheme, Geosci. Model Dev., 17, 1297–1326, https://doi.org/10.5194/gmd-17-1297-2024, 2024.

Bellaire, S., Katurji, M., Schulmann, T., and Hobman, A.: Towards a high-resolution operational forecasting tool for the southern alps – New Zealand, in: Proceedings, International Snow Science Workshop, 29 September 2014 to 3 October 2014, Banff, Alberta, Canada, https://arc.lib.montana.edu/snow-science/objects/ISSW14_paper_P1.01.pdf (last access: 24 June 2024), 2014.

Bergeron, T.: On the low-level redistribution of atmospheric water caused by orography, in: Proceedings of the International Conference on Cloud Physics, 24 May–1 June 1965, Tokyo and Sapporo, 96–100, https://doi.org/10.1088/0031-9112/16/9/004, 1965.

Bernhardt, M. and Schulz, K.: SnowSlide: A simple routine for calculating gravitational snow transport, Geophys. Res. Lett., 37, L11502, https://doi.org/10.1029/2010GL043086, 2010.

Bernhardt, M., Liston, G. E., Strasser, U., Zängl, G., and Schulz, K.: High resolution modelling of snow transport in complex terrain using downscaled MM5 wind fields, The Cryosphere, 4, 99–113, https://doi.org/10.5194/tc-4-99-2010, 2010.

Bernhardt, M., Schulz, K., Liston, G. E., and Zängl, G.: The influence of lateral snow redistribution processes on snow melt and sublimation in alpine regions, J. Hydrol., 424–425, 196–206, https://doi.org/10.1016/j.jhydrol.2012.01.001, 2012.

Brauchli, T., Trujillo, E., Huwald, H., and Lehning, M.: Influence of slope-scale snowmelt on catchment response simulated with the Alpine3D model, Water Resour. Res., 53, 10723–10739, https://doi.org/10.1002/2017WR021278, 2017.

Brun, E., Martin, E., Simon, V., Gendre, C., and Coleou, C.: An Energy and Mass Model of Snow Cover Suitable for Operational Avalanche Forecasting, J. Glaciol., 35, 333–342, https://doi.org/10.3189/S0022143000009254, 1989.

Candille, G. and Talagrand, O.: Evaluation of probabilistic prediction systems for a scalar variable, Q. J. Roy. Meteorol. Soc., 131, 2131–2150, https://doi.org/10.1256/qj.04.71, 2005.

Champeaux, J.-L., Dupuy, P., Laurantin, O., Soulan, I., Tabary, P., and Soubeyroux, J.-M.: Les mesures de précipitations et l'estimation des lames d'eau à Météo-France: état de l'art et perspectives, La Houille Blanche, 95, 28–34, https://doi.org/10.1051/lhb/2009052, 2009.

Choularton, T. W. and Perry, S. J.: A model of the orographic enhancement of snowfall by the seeder-feeder mechanism, Q. J. Roy. Meteorol. Soc., 112, 335–345, https://doi.org/10.1002/qj.49711247204, 1986.

Clark, M. P., Hendrikx, J., Slater, A. G., Kavetski, D., Anderson, B., Cullen, N. J., Kerr, T., Hreinsson, E. Ö., and Woods, R. A.: Representing spatial variability of snow water equivalent in hydrologic and land-surface models: A review, Water Resour. Res., 47, W07539, https://doi.org/10.1029/2011WR010745, 2011.

Cluzet, B., Lafaysse, M., Cosme, E., Albergel, C., Meunier, L.-F., and Dumont, M.: CrocO_v1.0: a particle filter to assimilate snowpack observations in a spatialised framework, Geosci. Model Dev., 14, 1595–1614, https://doi.org/10.5194/gmd-14-1595-2021, 2021.

CNRM/Centre d'Etudes de la Neige: The S2M meteorological and snow cover reanalysis in the French mountainous areas (1958–present), CNRM/Centre d'Etudes de la Neige [data set], https://doi.org/10.25326/37#v2020.2, 2019.

Colle, B. A., Smith, R. B., and Wesley, D. A.: Theory, observations, and predictions of orographic precipitation, in: Mountain Weather Research and Forecasting: Recent Progress and Current Challenges, Chap. 6, 291–344, Springer, https://doi.org/10.1007/978-94-007-4098-3_6, 2013.

Dadic, R., Mott, R., Lehning, M., and Burlando, P.: Wind influence on snow depth distribution and accumulation over glaciers, J. Geophys. Res.-Earth, 115, F01012, https://doi.org/10.1029/2009JF001261, 2010.
Decharme, B., Boone, A., Delire, C., and Noilhan, J.: Local evaluation of the Interaction between Soil Biosphere Atmosphere soil multilayer diffusion scheme using four pedotransfer functions, J. Geophys. Res., 116, D20126, https://doi.org/10.1029/2011JD016002, 2011.

Deschamps-Berger, C., Gascoin, S., Berthier, E., Deems, J., Gutmann, E., Dehecq, A., Shean, D., and Dumont, M.: Snow depth mapping from stereo satellite imagery in mountainous terrain: evaluation using airborne laser-scanning data, The Cryosphere, 14, 2925–2940, https://doi.org/10.5194/tc-14-2925-2020, 2020.

Deschamps-Berger, C., Cluzet, B., Dumont, M., Lafaysse, M., Berthier, E., Fanise, P., and Gascoin, S.: Improving the Spatial Distribution of Snow Cover Simulations by Assimilation of Satellite Stereoscopic Imagery, Water Resour. Res., 58, e2021WR030271, https://doi.org/10.1029/2021WR030271, 2022.

Durand, Y., Brun, E., Mérindol, L., Guyomarc’h, G., Lesafre, B., and Martin, E.: A meteorological estimation of relevant parameters for snow models, Ann. Glaciol., 18, 65–71, https://doi.org/10.3189/S0260305500011277, 1993.

Earth System Modeling Framework, University Corporation for Atmospheric Research, Massachusetts Institute of Technology, Geophysical Fluid Dynamics Laboratory, University of Michigan, National Centers for Environmental Prediction, Los Alamos National Laboratory, Argonne National Laboratory, and NASA Goddard Space Flight Center: Earth System Modeling Framework (ESMF) Grid Remapping (Version 8.4.2), https://earthsystemmodeling.org/ (last access: 24 June 2024), 2023.

Essery, R., Li, L., and Pomeroy, J.: A distributed model of blowing snow over complex terrain, Hydrol. Process., 13, 2423–2438, https://doi.org/10.1002/(SICI)1099-1085(199910)13:14/15<2423::AID-HYP853>3.0.CO;2-U, 1999.

Faure, D.: Quality Analysis of the 2016 Quantitative Precipitation Estimates in the French Alps, in: 37th Conference on radar meteorology, AMS, https://ams.confex.com/ams/38RADAR/webprogram/Paper320389.html (last access: 24 June 2024), 2017.

Faure, D., Delrieu, G., and Gaussiat, N.: Impact of the Altitudinal Gradients of Precipitation on the Radar QPE Bias in the French Alps, Atmosphere, 10, 306, https://doi.org/10.3390/atmos10060306, 2019.

Foresti, L., Sideris, I. V., Panziera, L., Nerini, D., and Germany, U.: A 10-year radar-based analysis of orographic precipitation growth and decay patterns over the Swiss Alpine region, Q. J. Roy. Meteorol. Soc., 144, 2277–2301, https://doi.org/10.1002/qj.3364, 2018.

Freudiger, D., Kohn, I., Seibert, J., Stahl, K., and Weiler, M.: Snow redistribution for the hydrological modeling of alpine catchments, WIREs Water, 4, e1232, https://doi.org/10.1002/wat2.1232, 2017.

Froidurot, S., Zin, I., Hingray, B., and Gautheron, A.: Sensitivity of precipitation phase over the Swiss Alps to different meteorological variables, J. Hydrometeorol., 15, 685–696, https://doi.org/10.1175/JHM-D-13-073.1, 2014.

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

Data, 11, 493–514, https://doi.org/10.5194/essd-11-493-2019, 2019.

GDAL/OGR contributors: GDAL/OGR Geospatial Data Abstraction software Library (Version 3.6.2), Open Source Geospatial Foundation, Zenodo [code], https://doi.org/10.5281/zenodo.5884351, 2023.

Germann, U., Boscacci, M., Clementi, L., Gabella, M., Hering, A., Sartori, M., Sideris, I. V., and Calpini, B.: Weather Radar in Complex Orography, Remote Sensing, 14, 503, https://doi.org/10.3390/rs14030503, 2022.

GLIMS and NSIDC: Global Land Ice Measurements from Space glacier database, International GLIMS community and the National Snow and Ice Data Center [data set], https://doi.org/10.7265/N5V98602, 2005.

Goessling, H. and Jung, T.: A probabilistic verification score for contours: Methodology and application to Arctic ice-edge forecasts, Q. J. Roy. Meteor. Soc., 144, 735–743, https://doi.org/10.1002/qj.3242, 2018.

Gouttevin, I., Vionnet, V., Seity, Y., Boone, A., Lafaysse, M., Deliot, Y., and Merzisen, H.: To the Origin of a Wintertime Screen-Level Temperature Bias at High Altitude in a Kilometric NWP Model, J. Hydrometeorol., 24, 53–71, https://doi.org/10.1175/JHM-D-21-0200.1, 2023.

Grünewald, T. and Lehning, M.: Are flat-field snow depth measurements representative? A comparison of selected index sites with areal snow depth measurements at the small catchment scale, Hydrol. Process., 29, 1717–1728, https://doi.org/10.1002/hyp.10295, 2015.

Grünewald, T., Stötter, J., Pomeroy, J. W., Dadic, R., Moreno Baños, I., Marturià, J., Spross, M., Hopkinson, C., Burlando, P., and Lehning, M.: Statistical modelling of the snow depth distribution in open alpine terrain, Hydrol. Earth Syst. Sci., 17, 3005–3021, https://doi.org/10.5194/hess-17-3005-2013, 2013.

Grünewald, T., Bühler, Y., and Lehning, M.: Elevation dependency of mountain snow depth, The Cryosphere, 8, 2381–2394, https://doi.org/10.5194/tc-8-2381-2014, 2014.

Günther, D., Marke, T., Essery, R., and Strasser, U.: Uncertainties in Snowpack Simulations – Assessing the Impact of Model Structure, Parameter Choice, and Forcing Data Error on Point-Scale Energy Balance Snow Model Performance, Water Resour. Res., 55, 2779–2800, https://doi.org/10.1029/2018WR023403, 2019.

Guyomarc’h, G. and Mérindol, L.: Validation of an application for forecasting blowing snow, Ann. Glaciol., 26, 138–143, https://doi.org/10.3189/1998AoG26-1-138-143, 1998.

Guyomarc’h, G., Bellot, H., Vionnet, V., Naaim-Bouvet, F., Déliot, Y., Fontaine, F., Puglièse, P., Nishimura, K., Durand, Y., and Naaim, M.: A meteorological and blowing snow data set (2000–2016) from a high-elevation alpine site (Col du Lac Blanc, France, 2720 m a.s.l.), Earth Syst. Sci. Data, 11, 57–69, https://doi.org/10.5194/essd-11-57-2019, 2019.

Haddjeri, A., Deschamp-Berger, C., and Gascoin, S.: Pleiades snow dataset from stereo-images, Zenodo [data set], https://doi.org/10.5281/zenodo.10037253, 2023.

Hofmeister, F., Arias-Rodriguez, L. F., Premier, V., Marin, C., Notarnicola, C., Disse, M., and Chiogna, G.: Intercomparison of Sentinel-2 and modelled snow cover maps in a high-elevation Alpine catchment, Journal of Hydrology X, 15, 100123, https://doi.org/10.1016/j.hydroa.2022.100123, 2022.
Horton, S. and Haegeli, P.: Using snow depth observations to provide insight into the quality of snowpack simulations for regional-scale avalanche forecasting, The Cryosphere, 16, 3393–3411, https://doi.org/10.5194/tc-16-3393-2022, 2022.

IGN©: BD Foret®V2, https://geoservices.ign.fr/bdforet, IGN© [data set], last access: 24 June 2024a.

IGN©: BD Topo®, https://geoservices.ign.fr/bdtopo, IGN© [data set], last access: 24 June 2024b.

IGN©: RGE ALTI®5M, https://geoservices.ign.fr/rgealti, IGN© [data set], last access: 24 June 2024c.

IGN©: BD Ortho®, IGN© [data set], https://geoservices.ign.fr/bdortho, 2022.

IPCC: The Ocean and Cryosphere in a Chang-
ing Climate, Cambridge University Press,
https://doi.org/10.1017/9781009157964, 2022.

Jasiewicz, J. and Stepinski, T. F.: Geomorphons – a pat-
tern recognition approach to classification and map-
ping of landforms, Geomorphology, 182, 147–156,
https://doi.org/10.1016/j.geomorph.2012.11.005, 2013.

Lafaysse, M., Cluzet, B., Dumont, M., Lejeune, Y., Vionnet, V., and Morin, S.: A multiphysical ensemble system of numerical snow modelling, The Cryosphere, 11, 1173–1198, https://doi.org/10.5194/tc-11-1173-2017, 2017.

Le Toumelin, L.: louisletoumelin/neural_network_and_devine: le_toumelin et al_2024 (le_toumelin_2024), Zenodo [code], https://doi.org/10.5281/zenodo.10594274, 2024.

Le Toumelin, L., Gouttevin, I., Helbig, N., Galiez, C., Roux, M., and Karbou, F.: Emulating the adaptation of wind fields to complex terrain with deep-learning, Artif. Intell. Earth Syst., 2, e220034, https://doi.org/10.1175/AIES-D-22-0034.1, 2022.

Le Toumelin, L., Gouttevin, I., Galiez, C., and Helbig, N.: A two-fold deep-learning strategy to correct and downscale winds over mountains, Nonlin. Processes Geophys., 31, 75–97, https://doi.org/10.5194/npg-31-75-2024, 2024.

Lehning, M., Völksch, I., Gustafsson, D., Nguyen, T. A., Stähli, M., and Zappa, M.: ALPINE3D: a detailed model of mountain surface processes and its application to snow hydrology, Hydrol. Process., 20, 2111–2128, https://doi.org/10.1002/hyp.6204, 2006.

Lehning, M., Löwe, H., Ryser, M., and Raderschall, N.: Inhomogeneous precipitation distribution and snow transport in steep terrain, Water Resour. Res., 44, W07404, https://doi.org/10.1029/2007WR006545, 2008.

Lindsay, J.: The whitebox geospatial analysis tools project and open-access GIS, in: Proceedings of the GIS Research UK 22nd Annual Conference, The University of Glasgow, 16–18, https://www.gla.ac.uk/media/Media_401757_smxx.pdf (last access: 24 June 2024), 2014.

Liston, G. and Sturm, M.: A snow-transport model for complex terrain, J. Glaciol., 44, 498–516, https://doi.org/10.3189/S0022143000002021, 1998.

Liston, G. E., Haehnel, R. B., Sturm, M., Hiemstra, C. A., Bere-
zovskaya, S., and Tabler, R. D.: Simulating complex snow distri-
butions in windy environments using SnowTran-3D, J. Glaciol.,
53, 241–256, https://doi.org/10.3189/172756507782202865,
2007.

Lundquist, J., Hughes, M., Gutmann, E., and Kapnick, S.: Our Skill in Modeling Mountain Rain and Snow is Bypassing the Skill of

Our Observational Networks, B. Am. Meteorol. Soc. 100, 2473–2490, https://doi.org/10.1175/BAMS-D-19-0001.1, 2019.

Marsh, C. B., Pomeroy, J. W., and Wheater, H. S.: The Canadian Hydrological Model (CHM) v1.0: a multi-scale, multi-extent, variable-complexity hydrological model – design and overview, Geosci. Model Dev., 13, 225–247, https://doi.org/10.5194/gmd-13-225-2020, 2020.

Marsh, C. B., Vionnet, V., and Pomeroy, J. W.: Windmap-
per: An Efficient Wind Downscaling Method for Hydro-
logical Models, Water Resour. Res., 59, e2022WR032683,
https://doi.org/10.1029/2022WR032683, 2023.

Marti, R., Gascoin, S., Berthier, E., de Pinel, M., Houet, T., and Laffly, D.: Mapping snow depth in open alpine terrain from stereo satellite imagery, The Cryosphere, 10, 1361–1380, https://doi.org/10.5194/tc-10-1361-2016, 2016.

Masson, V., Le Moigne, P., Martin, E., Faroux, S., Alias, A., Alkama, R., Belamari, S., Barbu, A., Boone, A., Bouyssel, F., Brousseau, P., Brun, E., Calvet, J.-C., Carrer, D., Decharme, B., Delire, C., Donier, S., Essaouini, K., Gibelin, A.-L., Giordani, H., Habets, F., Jidane, M., Kerdraon, G., Kourzeneva, E., Lafaysse, M., Lafont, S., Lebeaupin Brossier, C., Lemonsu, A., Mahfouf, J.-F., Marguinaud, P., Mokhtari, M., Morin, S., Pigeon, G., Salgado, R., Seity, Y., Taillefer, F., Tanguy, G., Tulet, P., Vincendon, B., Vionnet, V., and Voldoire, A.: The SURFEXv7.2 land and ocean surface platform for coupled or offline simulation of earth surface variables and fluxes, Geosci. Model Dev., 6, 929–960, https://doi.org/10.5194/gmd-6-929-2013, 2013.

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., and Wever, N.: 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, 865–880, https://doi.org/10.5194/essd-11-865-2019, 2019.

Menard, C. B., Essery, R., Krinner, G., Arduini, G., Bartlett, P., Boone, A., Brutel-Vuilmet, C., Burke, E., Cuntz, M., Dai, Y., Decharme, B., Dutra, E., Fang, X., Fierz, C., Gusev, Y., Hagemann, S., Haverd, V., Kim, H., Lafaysse, M., Marke, T., Na sonova, O., Nitta, T., Niwano, M., Pomeroy, J., Schädler, G., Se menov, V. A., Smirnova, T., Strasser, U., Swenson, S., Turkov, D., Wever, N., and Yuan, H.: Scientific and Human Errors in a Snow Model Intercomparison, B. Am. Meteorol. Soc., 102, E61–E79, https://doi.org/10.1175/BAMS-D-19-0329.1, 2021.

Minder, J. R., Durran, D. R., and Roe, G. H.: Mesoscale controls on the mountainside snow line, J. Atmos. Sci., 68, 2107–2127, https://doi.org/10.1175/JAS-D-10-05006.1, 2011.

Mital, U., Dwivedi, D., Brown, J. B., and Steefel, C. I.: Down-scaled hyper-resolution (400 m) gridded datasets of daily precipitation and temperature (2008–2019) for the East–Taylor subbasin (western United States), Earth Syst. Sci. Data, 14, 4949–4966, https://doi.org/10.5194/essd-14-4949-2022, 2022.

Monteiro, D., Caillaud, C., Samacoïts, R., Lafaysse, M., and Morin, S.: Potential and limitations of convection-permitting CNRM-AROME climate modelling in the French Alps, Int. J. Climatol., 42, 7162–7185, https://doi.org/10.1002/joc.7637, 2022.

Morin, S., Horton, S., Techel, F., Bavay, M., Coléou, C., Fierz, C., Gobiet, A., Hagenmuller, P., Lafaysse, M., Ližar, M., Miterer, C., Monti, F., Müller, K., Olefs, M., Snook, J. S., van Herwijnen, A., and Vionnet, V.: 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. Tech., 170, 102910,
https://doi.org/10.1016/j.coldregions.2019.102910, 2020.

Mott, R., Faure, F., Lehning, M., Löwe, H., Hynek, B., Michlmayer, G., Prokop, A., and Schöner, W.: Simulation of seasonal snow-cover distribution for glacierized sites on Sonnblick, Austria, with the Alpine3D model, Ann. Glaciol., 49, 155–160, https://doi.org/10.3189/17275640878714924, 2008.

Mott, R., Schirmer, M., Bavay, M., Grünewald, T., and Lehning, M.: Understanding snow-transport processes shaping the mountain snow-cover, The Cryosphere, 4, 545–559, https://doi.org/10.5194/tc-4-545-2010, 2010.

Mott, R., Scipión, D., Schneebeli, M., Dawes, N., Berne, A., and Lehning, M.: Orographic effects on snow deposition patterns in mountainous terrain, J. Geophys. Res.-Atmos., 119, 1419–1439, https://doi.org/10.1002/2013JD019880, 2014.

Mott, R., Vionnet, V., and Grünewald, T.: The Seasonal Snow Cover Dynamics: Review on Wind-Driven Coupling Processes, Front. Earth Sci., 6, 197, https://doi.org/10.3389/feart.2018.00197, 2018.

Mott, R., Winstral, A., Cluzet, B., Helbig, N., Magnusson, J., Mazzotti, G., Quéno, L., Schirmer, M., Webster, C., and Jonas, T.: Operational snow-hydrological modeling for Switzerland, Front. Earth Sci., 11, 1228158, https://doi.org/10.3389/feart.2023.1228158, 2023.

Musselman, K. N., Pomeroy, J. W., Essery, R. L., and Leroux, N.: Impact of windflow calculations on simulations of alpine snow accumulation, redistribution and ablation, Hydrol. Process., 29, 3983–3999, https://doi.org/10.1002/hyp.10595, 2015.

Pepin, N., Bradley, R. S., Diaz, H. F., Baraer, M., Caceres, E. B., Forsythe, N., Fowler, H., Greenwood, G., Hashmi, M. Z., Liu, X. D., Miller, J. R., Ning, L., Ohmura, A., Palazzi, E., Rangwala, I., Schöner, W., Severskiy, I., Shahgedanova, M., Wang, M. B., Williamson, S. N., Yang, D. Q., and Mountain Research Initiative EDW Working Group: Elevation-dependent warming in mountain regions of the world, Nat. Clim. Change, 5, 424–430, https://doi.org/10.1038/nclimate2563, 2015.

Pomeroy, J. and Gray, D.: Snow accumulation, relocation and management, vol. 144, National Hydrology Research Institute, https://hdl.handle.net/10388/15161 (last access: 24 June 2024), 1995.

Pomeroy, J., Gray, D., Shook, K., Toth, B., Essery, R., Pietroniro, A., and Hedstrom, N.: An evaluation of snow accumulation and ablation processes for land surface modelling, Hydrol. Process., 12, 2339–2367, https://doi.org/10.1002/(SICI)1099-1085(199812)12:15<2339::AID-HYP800>3.0.CO;2-L, 1998.

Pomeroy, J., MacDonald, M., DeBeer, C., and Brown, T.: Modelling Alpine Snow Hydrology in the Canadian Rocky Mountains, https://citeseer.xist.psu.edu/document?rid=rep1&type=pdf& doi=6ab63f9ac8702e66ff621b87c252ed6eaf025a6d (last access: 24 June 2024), 2009.

Prasad, R., Tarboton, D. G., Liston, G. E., Luce, C. H., and Seyfried, M. S.: Testing a blowing snow model against distributed snow measurements at Upper Sheep Creek, Idaho, United States of America, Water Resour. Res., 37, 1341–1356, https://doi.org/10.1029/2000WR900317, 2001.

Quéno, L., Vionnet, V., Dombrowski-Etchevers, I., Lafaysse, M., Dumont, M., and Karbou, F.: Snowpack modelling in the Pyrenees driven by kilometric-resolution meteorological forecasts,

The Cryosphere, 10, 1571–1589, https://doi.org/10.5194/tc-10-1571-2016, 2016.

Quéno, L., Karbou, F., Vionnet, V., and Dombrowski-Etchevers, I.: Satellite-derived products of solar and longwave irradiances used for snowpack modelling in mountainous terrain, Hydrol. Earth Syst. Sci., 24, 2083–2104, https://doi.org/10.5194/hess-24-2083-2020, 2020.

Quéno, L., Mott, R., Morin, P., Cluzet, B., Mazzotti, G., and Jonas, T.: Snow redistribution in an intermediate-complexity snow hydrology modelling framework, EGUSphere [preprint], https://doi.org/10.5194/egusphere-2023-2071, 2023.

Raleigh, M. S., Lundquist, J. D., and Clark, M. P.: Exploring the impact of forcing error characteristics on physically based snow simulations within a global sensitivity analysis framework, Hydrol. Earth Syst. Sci., 19, 3153–3179, https://doi.org/10.5194/hess-19-3153-2015, 2015.

Revuelto, J., Lecourt, G., Lafaysse, M., Zin, I., Charrois, L., Vionnet, V., Dumont, M., Rabatel, A., Six, D., Condom, T., et al.: Multi-Criteria Evaluation of Snowpack Simulations in Complex Alpine Terrain Using Satellite and In Situ Observations, Remote Sens., 10, 1171, https://doi.org/10.3390/rs10081171, 2018.

RGI Consortium: Randolph Glacier inventory-a dataset of global glacier outlines: Version 6.0, NSIDC [data set], https://doi.org/10.7265/4m1f-gd79, 2017.

Rouault, E., Warmerdam, F., Schwehr, K., Kiselev, A., Butler, H., Łoskot, M., Szekeres, T., Tourigny, E., Landa, M., Miara, I., Elliston, B., Chaitanya, K., Plesea, L., Morissette, D., Jolma, A., and Dawson, N.: GDAL (v3.7.2), Zenodo [code], https://doi.org/10.5281/zenodo.8340595, 2023.

Rüschendorf, L.: The Wasserstein distance and approximation theorems, Z. Wahrscheinlichkeitstheorie verw. Gebiete, 70, 117–129, https://doi.org/10.1007/BF00532240, 1985.

Sato, T., Kimura, T., Ishimaru, T., and Maruyama, T.: Field test of a new snow-particle counter (SPC) system, Ann. Glaciol., 18, 149–154, https://doi.org/10.3189/S0260305500011411, 1993.

Schirmer, M. and Jamieson, B.: Verification of analysed and forecasted winter precipitation in complex terrain, The Cryosphere, 9, 587–601, https://doi.org/10.5194/tc-9-587-2015, 2015.

Schlögl, S., Marty, C., Bavay, M., and Lehning, M.: Sensitivity of Alpine3D modeled snow cover to modifications in DEM resolution, station coverage and meteorological input quantities, Environ. Modell. Softw., 83, 387–396, https://doi.org/10.1016/j.envsoft.2016.02.017, 2016.

Schneiderbauer, S. and Prokop, A.: The atmospheric snow-transport model: SnowDrift3D, J. Glaciol., 57, 526–542, https://doi.org/10.3189/002214311796905677, 2011.

Seity, Y., Brousseau, P., Malardel, S., Hello, G., Bénard, P., Boutier, F., Lac, C., and Masson, V.: The AROME-France convective scale operational model, Mon. Weather Rev., 129, 976–991, https://doi.org/10.1175/2010MWR3425.1, 2011.

Sen Gupta, A. and Tarboton, D. G.: A tool for downscaling weather data from large-grid reanalysis products to finer spatial scales for distributed hydrological applications, Environ. Modell. Softw., 84, 50–69, https://doi.org/10.1016/j.envsoft.2016.06.014, 2016.

Sharma, V., Gerber, F., and Lehning, M.: Introducing CRYOWRF v1.0: multiscale atmospheric flow simulations with advanced snow cover modelling, Geosci. Model Dev., 16, 719–749, https://doi.org/10.5194/gmd-16-719-2023, 2023.
Strasser, U., Bernhardt, M., Weber, M., Liston, G. E., and Mauser, W.: Is snow sublimation important in the alpine water balance?, The Cryosphere, 2, 53–66, https://doi.org/10.5194/tc-2-53-2008, 2008.

Theia Snow: Theia Snow, Theia [data set], https://doi.org/10.24400/329360/F7Q52MNK, 2019.

Theurich, G., Collins, N., silveriovasquez10880, et al.: esmf-org/esmf: ESMF 8.6.1 (v8.6.1), Zenodo [software], https://doi.org/10.5281/zenodo.11205527, 2024.

Tong, J., Déry, S. J., and Jackson, P. L.: Topographic control of snow distribution in an alpine watershed of western Canada inferred from spatially-filtered MODIS snow products, Hydrol. Earth Syst. Sci., 13, 319–326, https://doi.org/10.5194/hess-13-319-2009, 2009.

Vernay, M., Lafaysse, M., Monteiro, D., Hagenmuller, P., Nheili, R., Samacoïts, R., Verfaillie, D., and Morin, S.: The S2M meteorological and snow cover reanalysis over the French mountainous areas: description and evaluation (1958–2021), Earth Syst. Sci. Data, 14, 1707–1733, https://doi.org/10.5194/essd-14-1707-2022, 2022.

Vionnet, V., Brun, E., Morin, S., Boone, A., Faroux, S., Le Moigne, P., Martin, E., and Willemet, J.-M.: The detailed snow-pack scheme Crocus and its implementation in SURFEX v7.2, Geosci. Model Dev., 5, 773–791, https://doi.org/10.5194/gmd-5-773-2012, 2012.

Vionnet, V., Martin, E., Masson, V., Guyomarc'h, G., Naaim-Bouvet, F., Prokop, A., Durand, Y., and Lac, C.: Simulation of wind-induced snow transport and sublimation in alpine terrain using a fully coupled snowpack/atmosphere model, The Cryosphere, 8, 395–415, https://doi.org/10.5194/tc-8-395-2014, 2014.

Vionnet, V., Dombrowski-Etchevers, I., Lafaysse, M., Quéno, L., Seity, Y., and Bazile, E.: Numerical weather forecasts at kilometer scale in the French Alps: evaluation and applications for snowpack modelling, J. Hydrometeorol., 17, 2591–2614, https://doi.org/10.1175/JHM-D-15-0241.1, 2016.

Vionnet, V., Martin, E., Masson, V., Lac, C., Naaim Bouvet, F., and Guyomarc'h, G.: High-Resolution Large Eddy Simulation of Snow Accumulation in Alpine Terrain, J. Geophys. Res.-Atmos., 122, 11005–11021, https://doi.org/10.1002/2017JD026947, 2017.

Vionnet, V., Six, D., Auger, L., Dumont, M., Lafaysse, M., Quéno, L., Réveillet, M., Dombrowski Etchevers, I., Thibert, E., and Vincent, C.: Sub-kilometer precipitation datasets for snowpack and glacier modeling in alpine terrain, Front. Earth Sci., 7, 182, https://doi.org/10.3389/feart.2019.00182, 2019.

Vionnet, V., Marsh, C. B., Menounos, B., Gascoin, S., Wayand, N. E., Shea, J., Mukherjee, K., and Pomeroy, J. W.: Multi-scale snowdrift-permitting modelling of mountain snowpack, The Cryosphere, 15, 743–769, https://doi.org/10.5194/tc-15-743-2021, 2021.

Vionnet, V., Verville, M., Fortin, V., Brugman, M., Abrahamowicz, M., Lemay, F., Thériault, J. M., Lafaysse, M., and Milbrandt, J. A.: Snow Level From Post-Processing of Atmospheric Model Improves Snowfall Estimate and Snowpack Prediction in Mountains, Water Resour. Res., 58, e2021WR031778, https://doi.org/10.1029/2021WR031778, 2022.

Wang, Z. and Huang, N.: Numerical simulation of the falling snow deposition over complex terrain, J. Geophys. Res.-Atmos., 122, 980–1000, https://doi.org/10.1002/2016JD025316, 2017.

Whitebox Geospatial Inc.: WhiteboxTools Open Core, Whitebox Geospatial Inc. [code], https://www.whiteboxgeo.com/manual/wbt_book/available_tools/geomorphometric_analysis.html# Geomorphons, last access: 24 June 2024.

Whiteman, C. D.: Mountain meteorology: fundamentals and applications, Oxford University Press, ISBN 9780197561522, https://doi.org/10.1093/oso/9780195132717.001.0001, 2000.

Winstral, A., Elder, K., and Davis, R. E.: Spatial snow modeling of wind-redistributed snow using terrain-based parameters, J. Hydrometeorol., 3, 524–538, https://doi.org/10.1175/1525-7541(2002)003<0524:SSMOWR>2.0.CO;2, 2002.

Yang, J., Yau, M. K., Fang, X., and Pomeroy, J. W.: A triple-moment blowing snow-atmospheric model and its application in computing the seasonal wintertime snow mass budget, Hydrol. Earth Syst. Sci., 14, 1063–1079, https://doi.org/10.5194/hess-14-1063-2010, 2010.

Yu, N., Gaussiat, N., and Tabary, P.: Polarimetric X-band weather radars for quantitative precipitation estimation in mountainous regions, Q. J. Roy. Meteor. Soc., 144, 2603–2619, https://doi.org/10.1002/qj.3366, 2018.