Canopy structure, topography, and weather are equally important drivers of small-scale snow cover dynamics in sub-alpine forests

. In mountain regions, forests that overlap with seasonal snow mostly reside in complex terrain. Due to persisting major observational challenges in these environments, the combined impact of forest structure and topography on seasonal snow cover dynamics is still poorly understood. Recent advances in forest snow process representation and increasing availability of detailed canopy structure datasets, however, now allow for hyper-resolution ( < 5 m) snow model simulations capable of resolving tree-scale processes. These can shed light on the complex process interactions that govern forest snow dynamics. We present multi-year simulations at 2 m resolution obtained with FSM2, a mass-and energy-balance-based forest snow model speciﬁcally developed and validated for metre-scale applications. We simulate an ∼ 3 km 2 model domain encompassing forested slopes of a sub-alpine valley in the eastern Swiss Alps and six snow seasons. Simulations thus span a wide range of canopy structures, terrain characteristics, and meteorological conditions. We analyse spatial and temporal variations in forest snow energy balance partitioning, aiming to quantify and understand the contribution of individual energy exchange processes at different locations and times. Our results suggest that snow cover evolution is equally affected by canopy structure, terrain characteristics, and meteorological conditions. We show that the interaction of these three factors can lead to snow accumulation and ablation patterns that vary be-tween years. We further identify higher snow distribution variability and complexity in slopes that receive solar radiation early in winter. Our process-level insights corroborate and complement existing empirical ﬁndings that are largely based on snow distribution datasets only. Hyper-resolution simulations as presented here thus help to better understand how snowpacks and ecohydrological regimes in sub-alpine regions may evolve due to forest disturbances and a warming climate. They could further support the development of process-based sub-grid forest snow cover parameterizations or tiling approaches for coarse-resolution modelling applications.


S1: Physiographic maps over the full model domain
To complement the site description in the main article (Section 2.2 and Figure 1), Figure S1.1 reports additional maps of physiographic variables across the model domain, including local canopy cover fraction, elevation, aspect and slope.The diversified representation of canopy structure of FSM2 with process-specific canopy metrics constitutes the model's principal asset.With this strategy, FSM2 considers that different canopy-mediated processes are regulated by different and potentially uncorrelated canopy characteristics.A point located in a canopy gap, as example, experiences the interception regime of a very sparse canopy, but may be affected by frequent shading, which implies conditions typical of a dense canopy.These independent and contrasting canopy controls on each of the various processes cannot be reproduced if all processes are parametrized with the same, bulk canopy descriptors.Bulk canopy metrics may thus be suited for intermediateor coarse-resolution model applications targeting the estimation of spatially averaged fluxes and states, but hyper-resolution approaches that strive to reproduce meter-scale patterns require explicit representation of a wide range of detailed canopy features.The canopy metrics used within FSM2 comprise horizontal, vertical, local, and stand-scale canopy characteristics at each modelled location.Consequently, the canopy's structural diversity is captured even though canopy is represented with a so-called one-layer model.Figure S2.1 illustrates this concept, summarizing the different canopy perspectives, the corresponding metrics, approaches, and data sources used to derive them, and the processes employing them.The schematic complements the mathematical formulation of the process parametrizations provided in the respective model description papers, Essery (2015) and Mazzotti et al. (2020a,b).
Transmission of direct shortwave radiation through the canopy involves a particularly complex interaction with small-scale canopy elements.Because transmissivity is dictated by the presence of canopy elements in the path of the solar beam, it depends on the exact geometric arrangement of the point of interest, the canopy, and the sun, and is thus highly variable in space and time.By accepting transmissivity time series as model input, FSM2 forgoes simplification of the process representation that comes with any parametrization, and instead leverages the benefits of an (external) radiative transfer model that resolves the process explicitly.So far, FSM2.0.3 has been used in tandem with the radiative transfer model HPEval (Jonas et al., 2020), which comprises a traditional approach to compute transmissivity time series at a point based on a high-resolution hemispherical image.The radiative transfer model depicted in the schematic in Figure S2.1 illustrates the HPEval approach.In this study, transmissivity time series are computed following the HPEval methodology based on synthetic hemispherical images derived from LiDAR datasets instead of real ones.The workflow to create these synthetic images, detailed in Webster et al. (2020), integrates canopy segmentation algorithms to identify individual trees in a canopy height model and point cloud enhancing algorithms to densify the point cloud, aimed at mimicking opaque stems and branches.The resulting images feature more realistic and detailed tree shapes than images obtained with unprocessed LiDAR data (see Section S3).

S3: Summary of FSM2 validation at the process level
This study relies on FSM2's capability to accurately simulate the spatio-temporal patterns of individual fluxes and states that constitute forest snow processes, which was extensively tested in a previous study by Mazzotti et al. (2020a).Figure S3.1 reports an example of their validation efforts to illustrate the underlying methodology.Data were acquired with multiple meteorological sensors mounted on a motorized cable car platform (Figure S3.1).Over the course of a winter, the system was deployed along several forest transects of approx.50 m length each in both Switzerland and Finland.Surveys covered variable canopy structures and forest species, weather conditions, as well as sub-alpine and boreal climates.Observations from the cable car setup were complemented by data from two additional handheld and stationary setups, yielding a dataset with an unprecedented variety of micro-meteorological and snow measurements.FSM2 simulations at high spatial and temporal resolution were directly compared to the acquired datasets.The example in Figure S3.1 features measurements from a full 24-hour cycle along a cable car transect and corresponding FSM2 simulations at 1.5m spacing and 10min temporal resolution.To our knowledge, FSM2 is the only forest snow model to have ever undergone such a rigorous validation.
FSM2's ability to replicate individual fluxes and states and their spatial and temporal variability was quantitatively assessed in terms of transect averages and standard deviations of all measured variables across all available datasets (Figure S3.2).
Root means square errors (RMSE) of 23 Wm −2 and 4Wm −2 resulted for average shortwave and longwave irradiances, respectively, while their standard deviations featured RMSEs of 21Wm −2 and 2Wm −2 .RMSEs of mean and standard deviation of canopy air space temperatures were 1.6°C and 0.4°C, for snow surface temperatures they amounted to 0.6°C and 0.4°C.Performance metrics all showed substantial improvements relative to the standard land surface model canopy implementation that was used as benchmark in Mazzotti et al. (2020a).
The accuracy of the shortwave radiative transfer calculations in FSM2 is dictated by the accuracy of the external radiative transfer model.Webster et al. (2020) showed that shortwave radiation estimates using synthetic images obtained with the LiDAR enhancing methodology were similar to estimates obtained with real images (Figure S3.3, a-e), and patterns of shadows on the snow surface (f) matched those captured with UAV-imagery (g) with high level of detail.
While plausibility checks as presented in Section 2.3 (main article) are important to ensure model applicability for a specific use case, we believe that quantitative validations at the level of individual fluxes and states such as in Mazzotti et al. (2020a) and summarized here are more conclusive than assessments based on snow extent and depth data alone.

S4: Meteorological conditions during the simulation period
The following figures provide an overview of meteorological conditions during the six simulated water years to give context to the observed differences in snow cover dynamics outlined in the main article.We focus on the three variables that, as elaborated on in the main article, have the largest potential to cause differences between years: snowfall (affects total accumulation), shortwave radiation, and temperature (affect the timing and strength of melt, as well as whether it is mainly driven by short-or longwave radiation contributions).The left panel of Figure S4.1 reports weekly snowfall sums between the beginning of November and June, which represent the earliest and latest yearly median of the start of snow cover period and snow disappearance date.Some marked differences in how the precipitation is distributed throughout the snow season can be observed.The right panel of Figure S4.1 shows cumulative snowfall between the median start of snow cover and snow disappearance dates of each specific year.This visualization evidences the substantial differences in timing and length of snow cover period and total snowfall between the years, with WY2016 (lowest accumulation) amounting to only half of total snowfall of WY 2019 (highest accumulation).

Figure S2. 1 :
Figure S2.1:Overview of the concepts underlying the canopy structure representation in FSM2: different perspectives, relevant data sources and computation approaches, resulting metrics, and modelled processes that use the respective metrics.

Figure S3. 1 :
Figure S3.1:The multi-sensor platform in operation along a 50-m cable (yellow arrow) in a pine stand near Davos, Switzerland 90

Figure S4. 1 :
Figure S4.1:Weekly snowfall during the six modelled winters between November and June (left) and cumulative snowfall over the median snow cover duration period of the respective year (right).

Figure
Figure S4.2 contrasts shortwave radiation (left) and air temperatures (right) during the six simulated winters, with the upper row showing weekly values and the lower row the weekly deviation from the six-year average.While no year is clearly 125above or below average for any of the two variables for the entire season, some marked deviations from the averages can be detected in these plots, for instance the relatively warm period in the middle of the 2020 winter.

Figure S4. 2 :
Figure S4.2:Incoming shortwave radiation (left) and air temperature (top right) of the six simulated winters, and deviation from mean (lower row) for the same variables.130

Figure
Figure S5.2 shows average incoming irradiances and early melt between mid-January and end of February 2019, analogous 140 to Figure 7 of the main article.Note the addition of incoming longwave radiation which further underpins the statement that all-wave irradiance patterns are largely determined by patterns of incoming shortwave radiation.

Figure S5. 2 :
Figure S5.2:Average incoming short (a) and longwave (b) radiation, all-wave radiation (c) and snow melt (d) between mid-January and end of February 2019.145