Snow farming : conserving snow over the summer season

Summer storage of snow for tourism has seen an increasing interest in the last years. Covering large snow piles with materials such as sawdust enables more than two-thirds of the initial snow volume to be conserved. We present detailed mass balance measurements of two sawdust-covered snow piles obtained by terrestrial laser scanning during summer 2015. Results indicate that 74 and 63 % of the snow volume remained over the summer for piles in Davos, Switzerland and Martell, Italy. If snow mass is considered instead of volume, the values increase to 83 and 72 %. The difference is attributed to settling and densification of the snow. Additionally, we adapted the one-dimensional, physically based snow cover model SNOWPACK to perform simulations of the sawdust-covered snow piles. Model results and measurements agreed extremely well at the point scale. Moreover, we analysed the contribution of the different terms of the surface energy balance to snow ablation for a pile covered with a 40 cm thick sawdust layer and a pile without insulation. Short-wave radiation was the dominant source of energy for both scenarios, but the moist sawdust caused strong cooling by long-wave emission and negative sensible and latent heat fluxes. This cooling effect reduces the energy available for melt by up to a factor of 12. As a result only 9 % of the net short-wave energy remained available for melt. Finally, sensitivity studies of the parameters “thickness of the sawdust layer”, “air temperature”, “precipitation” and “wind speed” were performed. We show that sawdust thickness has a tremendous effect on snow loss. Higher air temperatures and wind speeds increase snow ablation but less significantly. No significant effect of additional precipitation could be found as the sawdust remained wet during the entire summer with the measured quantity of rain. Setting precipitation amounts to zero, however, strongly increased melt. Overall, the 40 cm sawdust provides sufficient protection for midelevation (approx. 1500 m a.s.l.) Alpine climates and can be managed with reasonable effort.


Introduction
Snow storage or snow farming is the conservation of snow during the warm season of the year.The purposes of summering snow or ice are diverse.In ancient times, ice and snow were stored to cool food and houses (Morley, 1942;Taylor, 1985;Skogsberg and Nordell, 2001;Morofsky, 2007).Another example of a traditional application of snow storage is the collection of snow in deep underground wells, e.g. in Afghanistan (Bhattacharyya et al., 2004).The water stored in the snow could be used for irrigation or as drinking water during summer.Rising energy costs have increased interest in snow or ice as a cooling source (Nordell, 2014).In several regions of the world, such as Scandinavia, North America, China or Japan snow storage has seen a revival for air conditioning of buildings (Skogsberg and Nordell, 2001;Nordell and Skogsberg, 2007;Morofsky, 2007;Hamada et al., 2010Hamada et al., , 2012) ) and food cooling (Kobiyama et al., 1997;Nordell, 2014): during winter, large amounts of natural or machinemade snow are collected and stored.During summer the melting snow provides cooling energy for air conditioning.
Another application of snow conservation is the protection of glaciers from melt.At some ski resorts, thin textile covers called geotextiles have been used to blanket critical areas of the glacier surface at the end of the winter season (Olefs and Lehning, 2010).The properties of the covering material, such as albedo, thermal conductivity, emissivity or permeability influence the energy balance of the snow and the glacier, resulting in reduced ablation (Olefs and Obleitner, 2007;Olefs and Lehning, 2010).Olefs and Fischer (2008) performed field tests at two Austrian glaciers where they analysed different types of geotextiles and plastic fabrics as covering material.They found that ablation could be reduced by up to 60 % in comparison to a natural snow surface.Finally, Olefs and Lehning (2010) used the snow cover model SNOWPACK to simulate the effect of geotextiles on energy balance and snow ablation on glaciers.They showed that most of the performance of geotextiles is attributed to increased short-wave reflectance compared to the dirty firn surface.Thermal insulation and cooling by latent heat transfer also contribute but less significantly.Olefs and Lehning (2010) summarize that geotextiles are very effective covering materials for glaciers but are not suitable at lower elevations where turbulent fluxes typically dominate the energy balance.
The idea to conserve snow for touristic purposes, also at lower elevations, came up in Scandinavia more than a decade ago.Similarly to the previously mentioned snow storage applications, large amounts of snow are collected at the end of the winter or produced by snow machines and conserved over the summer months.In the following autumn, the remaining snow is used as basis for the preparation of winter sport facilities such as cross-country tracks, ski runs or ski jumps.A frequent motivation of winter sport destinations to store snow is the hosting of large sports events, such as crosscountry races in autumn or early winter.The conserved snow provides a guarantee for a basic amount of snow for track preparation, independent of the current weather conditions, which may be unfavourable for technical snow production.Examples of such events are the Cross-Country World Cup in Davos, Switzerland, the Biathlon World Cup in Östersund, Sweden or the Ski Jumping World Cup in Titisee-Neustadt, Germany.The largest application of snow storage was performed for the winter Olympics held in Sotschi, Russia in 2014.About 800 000 m 3 of snow was preserved as a reserve for the preparation of alpine ski-racing tracks in case of a lack of snow (Lintzen, 2016).Other motivations are the early preparation of training facilities for athletes or an early opening of ski slopes for the public.Some ski resorts collect snow at the end of the winter and conserve it at neuralgic spots, such as terrain breaks or depressions.In the succeeding winter the remaining snow is used to level these depressions.
In recent years the number of snow conservation projects and sites has strongly increased, especially in Scandinavia but also in the Alps.Snow farming can be seen as an adaption to climate warming but is also required to satisfy consumer demands of an early winter season start.However, the number of additional tourists in the early winter season due to snow farming is probably limited (Dreier, 2010).The main benefit for the destinations will rather be related to a higher media presence and a positive image in relation to general snow security of the destination and for hosting professional winter sport events such as World Cup races regardless of the weather conditions.
Different organic materials, such as sawdust, chipped wood, cutter shaving, bark mulch or straw have been used to cover the snow piles (Skogsberg and Lundberg, 2005).Moreover, non-organic fabrics such as styrofoam plates, polymeric foam, geotextiles and combinations of different materials have been applied.All these materials act as insulating layers that reduce heat transfer from the atmosphere to the snow.Depth, thermal conductivity and heat capacity are the key characteristics for their insulation performance.In addition, most organic materials are able to store significant amounts of moisture.Evaporation of the water cools the surface and reduces snow ablation (Skogsberg and Lundberg, 2005).Moreover, the high short-wave reflectivity (albedo) of some cover materials reduces solar energy input.Obviously, meteorological conditions, such as air temperature, solar radiation, humidity, precipitation and wind strongly affect snowmelt.First investigations (Rinderer, 2009), basic calculations (Skogsberg and Nordell, 2001) and reports from practitioners indicate that the amount of snow that is lost during summer is in the range of 20 to 50 %.However, mass balances based on accurate measurements have not yet been published.The only study known to the authors is a field test by Rinderer (2009) who observed two snow heaps at a location near Davos, CH (1620 m a.s.l.).One was covered with geotextile (4 mm) and the other with a 40 cm layer of sawdust.Based on a simple analysis of time lapse photography, Rinderer (2009) concluded that the first pile nearly melted completely while the latter lost about 30 % of its volume.Additionally, the different cover layers were implemented in the physically based snow cover model SNOWPACK (Lehning et al., 2002a).Simulations of the field experiments showed plausible results, but due to the limited quantification of the volume losses a precise validation of the model results were lacking.
A rising number of technical feasibility studies on snow farming has been requested at SLF in recent years.This motivated us to (i) provide a review on current snow-farming techniques (ii) to perform a detailed field study on snow farming and (iii) to describe and evaluate the model used for snow-farming applications by the SLF.This paper presents the first scientific publication on snow farming for touristic purposes as described above.Contrary to Rinderer's relatively rough estimation, we present detailed measurements of the volume loss based on high-resolution terrestrial laser scanning surveys.In 2015, measurement campaigns were performed at two sites, at the same site as Rinderer's study (Davos) and in the Martell valley of South Tyrol.We present an analysis of the small-scale spatial characteristics of snow ablation at the two snow heaps.The study is completed by simulations using the SNOWPACK model.The model was run for the two sites to evaluate the ability of SNOWPACK to reproduce the evolution of the snow piles and to analyse the impact of the different terms of the surface energy balance.
Finally, a sensitivity study, assessing effects of different settings and changed atmospheric conditions is performed.The respective results are discussed and summarized at the end of the paper.
2 Methods and data

Study sites and snow deposits
Observations of snow-farming projects at two study sites in the Alps are presented.The first is in the Flüela valley near the city of Davos in the eastern part of Switzerland (Fig. 1).The storage site is a forest clearing at the valley bottom located at 1620 m a.s.l.near to the mouth of the Flüela valley (lat: 46.808 • N, long: 9.868 • E).The direction of the valley is west to east and the site is shadowed by mountains, in particular by a 2500 m high mountain in the South.
Snow farming has been operated at the Flüela site since 2008.A large snow pile is formed by machine-made snow produced during the winter months.The snow pile also contains relatively small amounts of natural snow resulting from precipitation events during winter.In spring this pile is smoothed and covered by a 40 cm layer of sawdust.In 2015 the snow heap was covered between 13 and 18 April and the sawdust was removed on 19 October.The purpose of snow farming in Davos is to provide a short track for cross-country skiing, which is usually open by the beginning of November.The track serves as training facility for professional athletes and also provides the snow basis for the Cross-Country World Cup in early winter.
Meteorological data (Table 1) from two different automatic stations were analysed for the Flüela site.The Windtunnel station, measuring air temperature (TA), relative humidity (RH), wind speed (WS) and direction (WD), is installed on the roof of a concrete building next to the snow heap.The height of the sensors is about 7 m above ground, which is similar to the height of the snow pile.The remaining meteorological parameters required as input for the simulations (incoming short-wave radiation, ISWR; incoming longwave radiation, ILWR; and precipitation, P ) were taken from the weather station DAV, which is located at a similar elevation (1596 m a.s.l.) in the town of Davos about 2 km away and has high-quality meteorological observations.The local effect of terrain shadowing to incoming short-wave radiation was corrected by applying a shade-filter available in the preprocessing library for meteorological data, MeteoIO (Bavay and Egger, 2014).This filter transfers measured radiation of a measurement station to another site taking into account local shading.
No automatic station had been set up directly on top of the snow pile in 2015.However, a simple weather station had been placed at the top of the snow heap from 27 May to 17 July 2016.These observations could be applied to assess the representativeness of the sensors and to correct the  The second snow pile was in the Martell valley of South Tyrol (Fig. 2).The site is a forest clearing at 1710 m a.s.l. at the bottom of the valley (lat: 46.517 • N, long: 10.727 • E).Up to 3700 m high mountains surround a valley that opens to the north-east.The snow-dumping site is a concave, northerly exposed slope, steepening from 0 to 50 • .The work flow and purpose of the pile are similar to Flüela: snow is produced with snow machines during the cold season and covered in spring.In contrast to Flüela, the covering material at Martell also contained a portion of larger wood chips.In 2015 the snow heap was covered from 6 May to 13 November: 3.5 km of cross-country tracks could be prepared with the stored snow in that year.Similarly to Davos, the aims of snow farming are the provision of early training possibilities for athletes and a guarantee of basic snow amounts, independent of weather conditions during autumn and early winter.
We used meteorological observations from the automatic weather station Hintermartell (Autonome Provinz Bozen, 2016b) to assess meteorological forcing.This station is only www.the-cryosphere.net/12/385/2018/The Cryosphere, 12, 385-400, 2018 250 m away and is equipped with sensors for TA, RH, P , ISWR, WS and WD.Measured ILWR is not available in the region and was therefore parametrized (see Sect. 2.4).
As the station Hintermartell was only established in 2009, we employed data from the slightly higher (1851 m a.s.l.) and 1.2 km distant station, Stausee Zufritt (SZ) (Autonome Provinz Bozen, 2016a) to evaluate the climate conditions of summer 2015 in a climatological context.

Terrestrial laser scanning
Measurements of the snow volume stored at the end of the winter season and of the volume that remained in autumn were performed by terrestrial laser scanning (TLS).TLS proved as a highly accurate method to obtain digital surface models.Changes in surface volume can easily be obtained by subtracting two succeeding surveys.In past studies repeated TLS was successfully applied to calculate snow volumes at high spatial resolution and with a vertical accuracy of less than 10 cm (e.g.Prokop et al., 2008;Grünewald et al., 2010;Revuelto et al., 2014).
We used a Riegl VZ6000 terrestrial laser scanner (Riegl Laser Measurement Systems GmbH, 2015) shown in Fig. 2. The VZ6000 operates in the near-infrared spectral range (1064 nm) and is characterized by high accuracy (15 mm), precision (10 mm) and a beam divergence of only 0.12 mrad.Measurement frequencies of 300 kHz were used.For this study, the maximum measurement distance of the target area was less than 150 m.The respective beam width for this distance would be 33 mm.Because of this rather small measurement range the angular step width could be reduced to 0.03 • .To overcome scan shadows and to capture the com-plete snow piles, TLS measurements from six (five) positions were taken at each survey day at the Flüela (Martell) snow pile.TLS surveys were performed at the Flüela site on 29 April and 8 October 2015.These two data sets allow for the calculation of the volume change of the snow heap during summer.In order to determine the absolute volume, and through this the relative loss of the snow pile, a bare-ground (without snow pile) elevation model would be required.This model was not available and could not be surveyed as the site was always covered, either by the snow heap (during summer) or by the stored sawdust (after snow removal).However, the site is a quite homogeneous, only slightly sloped flat area that could well be approximated by a plane interpolated from the edges of the snow pile.Note that at the time of both surveys the snow pile was covered by a layer of sawdust.In autumn 23 manual measurements of the sawdust depth were performed in two transects with an avalanche probe perpendicular to the slope.The mean depth was 32 cm (standard deviation 11 cm), corresponding to a depth of about 38 cm when calculated in the vertical (gravitational) direction.
At the Martell field site TLS measurements were performed on 19 May and 28 October 2015.Similarly to Flüela these dates indicate the beginning (a few days after the pile was covered) and the end of the melt season (briefly before the snow was turned out to the cross country track).Identical to the Flüela data set the period between the two surveys corresponds to 163 days but both spring and autumn measurements were obtained 3 weeks later at Martell.In contrast to Flüela, the covering material was not only sawdust but sawdust with wood chips.The particle sizes of these wood chips varied from 0.1 mm (like sawdust) to several centimetres.It needs to be considered that the surface of the pile was additionally covered by a thin layer of snow during the second survey that resulted from an early snow fall.We obtained a mean snow height (HS) of 6 cm (standard deviation 1.5 cm) from 27 manual HS measurements.During data processing, this mean value was applied to remove the additional volume caused by the snow cover.As the wood chips layer was hard frozen, its thickness could not be measured with a penetrating probe but, based on the roughly known total volume of wood chips (800-1000 m 3 ) deployed it was estimated to be about 35 to 45 cm (Martin Stricker personal communication).As mentioned before, the ground surface of the dump is not flat at Martell but characterized by sloped topography.Therefore, a third TLS survey was performed on 19 November.Even though most snow had been turned out at that time, some snow still remained in the dump.Unfortunately, HS could not be measured with a probe as the snow had frozen hard.We therefore had to estimate the remaining snow amount based on our visual impression: we subdivided the area into three zones, one with estimated HS of 10 to 20 cm, one with 50-100 cm and one with 150-200 cm.These zones were then mapped with a differential GPS.A potential range of snow volume remaining in the dump could be estimated (600 to 1000 m 3 ).Finally, we used a mean value of 800 m 3 for the corrections of the heap volumes.

Snow density
Snow densities are, on the one hand, required for the initialization of the SNOWPACK model and, on the other hand, for the transformation from snow volume to snow mass or snow water equivalent (SWE).Snow densities were measured at different locations and depths of the Flüela snow heap at the end of the storage period (October 2015) and at the snow heap of the succeeding year (April 2016).In spring, five snow volumes were extracted with a core driller providing cylindric snow samples (d = 72 mm) of different length.In autumn, 10 volume samples were taken with a 100 dm 3 cubic density cutter.Additionally, three cubic snow samples were taken to determine densities by X-ray micro-computer tomography (micro-CT) (Heggli et al., 2011).No densities could be measured at Martell.

Data processing
TLS required substantial efforts of data postprocessing, which was done with the commercial software RiSCAN Pro provided by the manufacturer (Riegl Laser Measurement Systems GmbH, 2011).In a first step the large amount of data needed to be reduced.This was done by removing all data outside of the area of interest and by aggregating the remaining data points to a 5 cm 3-D grid (octree filter).Then the data points from the multiple-scan positions needed to be combined.For this registration, a set of eight (seven) reflector tiepoints were installed at fixed locations such as the walls of buildings in the surrounding.Global coordinates of the reflectors were recorded with differential GPS and total station.In a second step an optimization function called multi-station adjustment was applied to further improve the registration of the scan positions (Riegl Laser Measurement Systems GmbH, 2011).Finally, data from each survey were transformed to the global coordinate system (Flüela: CH1903; Martell: UTM) by applying the global coordinates of the reflectors and the merged and filtered point clouds were exported.ESRI ArcMap 10.2 was used to rasterize the data.Data were triangulated and data gaps that were partly existing, especially at the crown of the pile, were closed by triangulation with the nearest points.Finally, raster data with a resolution of 10 cm were composed and used for volume calculations and further analysis.

Snowpack model
SNOWPACK is a one-dimensional physically based snow model which has been developed to simulate state and evolution of the snow cover.The multilayer model includes a detailed description of the energy and mass fluxes within the snow and between snow, atmosphere and soil (e.g.Lehning et al., 2002a, b;Wever et al., 2015).For this publica-tion SNOWPACK was applied to simulate the evolution of the two snow heaps during a complete summer season.The aims were to identify the dominant processes related to mass conservation during snow farming and to evaluate the models applicability to improve and plan existing and upcoming snow-farming projects.
SNOWPACK is initialized with a snow profile, describing the state of the snow at the beginning of the simulation.With the one-dimensional SNOWPACK model, we chose to simulate ablation at the top of the respective piles.Corresponding to the maximal HS measured at the two heaps, HS was set to 8.60 m for Flüela and 7.20 m for Martell.The snow was subdivided into 10 homogeneous layers of fully rounded grains with a grain radius of 1 mm and a density of 553 kg m −3 .These values are typical for settled machine-made snow and are based on the measurements described earlier.Moreover, the snow was assumed to be isothermal (snow temperature 0 • C), holding a portion of 3 % liquid water.Temperature and water content were not measured directly but are based on our impressions during density sampling.Two extra layers were added to the snow profile, one at the bottom (ground) and one at the surface representing the covering material.A 40 cm thick covering layer representing the amount of sawdust obtained from measurements at the two snow heaps is used as reference simulation in the following.Characteristics of the sawdust layer are listed in Table 2. SNOWPACK allows the representation of soil (Lütschg et al., 2008) and other materials (e.g.Olefs and Lehning, 2010) by specifying mechanical and thermal layer properties.As we assume similar properties of sawdust and the mixture of sawdust and wood chips used in Martell, the same model settings were used for all simulations.
SNOWPACK is driven by the meteorological measurements TA, RH, WS, WD, ISWR, ILWR and P (Table 1).These parameters were obtained from nearby automatic weather stations as described earlier.As no measurements for ILWR were available for the Martell region, it was approximated using a combination of a clear-sky parametrization (Dilley and O'Brien, 1998) and an all-sky algorithm (Unsworth and Monteith, 1975) as described in Schmucki et al. (2014).This parametrization was also applied to account for changes in long-wave radiation that go with increasing temperatures (Schmucki et al., 2014) as simulated in the sensitivity runs with higher air temperatures (Sect.3.5).
All input data were filtered, quality checked and resampled to the modelling time step of 15 min using the meteorological input-output library MeteoIO (Bavay and Egger, 2014).Model outputs are time series of snow profiles and fluxes reflecting the state of the snowpack at different points in time.In the context of this study, snow height, snow mass, density, water content and the respective energy and mass fluxes are the parameters of interest.In addition to the reference simulation with measured meteorological data and a 40 cm sawdust cover, sensitivity studies were performed by varying the meteorological input parameters and the depth www.the-cryosphere.net/12/385/2018/The Cryosphere, 12, 385-400, 2018 of the covering surface layer (Table 3).For each run, only one of the parameters was changed while the measured input was used for the others (as described above).We performed seven sensitivity runs for different depths of the covering material by increasing it from 0 to 60 cm.In addition to the reference simulation, we performed two model runs for wind speed, three runs for precipitation and four runs for temperature as shown in Table 3.Note that ILWR was also adapted for the temperature simulations as described above, which increases the consistency of the temperature sensitivity runs (Schmucki et al., 2014).Sensitivity runs are, on the one hand, valuable for identifying the most relevant impact factors that can help to improve settings for snow-farming projects and, on the other hand, test the operability of snow farming as presented to other sites or under future climatic conditions.

Meteorological data
Half-hourly meteorological measurements of WS and TA obtained by the stations at the top of the Flüela snow pile and the station Windtunnel (WT) were analysed for the period 27 May to 17 July 2016.Air temperature at Windtunnel was highly correlated (R 2 = 0.9) and showed a low bias of 0.4 • C with measurements at the pile.Contrary to TA, WS showed stronger deviations.A mean WS of 0.3 m s −1 measured at Windtunnel indicates a clear underestimation of the wind at the pile (mean = 0.9 m s −1 ).This might be attributed to wind-sheltering effects of local topography and surrounding trees.Based on the measured data we calculated a linear regression function (WS = 2.1 × WS WT + 0.3, R 2 = 0.6) that was then applied to correct WS for the simulations.
Meteorological input data applied to the simulations are presented in Figs. 3 and 4. Note that the data were aggregated to daily values for the plot while half-hourly values were used as input for the simulations.Mean temperatures of the survey period (29 April to 8 October) were 11.Precipitation cumulated to 571 mm (SZ: 409 mm) in comparison to a 15-year mean of 645 mm (SZ: 434 mm).It can therefore be concluded that the summer half of the year 2015 was warmer and drier at both sites, while the difference was more pronounced at Flüela.This is in agreement with the findings of climate reports for Switzerland (MeteoSchweiz, 2016) and South Tyrol (Munarni et al., 2015) that rated the year 2015 as the warmest and the summer as the second warmest since the start of the measurements.

Snow density
Four of the five density samples collected at the snow pile in spring could be used for micro-CT analysis.One sample had been destroyed during preparation.Mean density calculated from the samples with a volume of about 5 cm 3 was 555 kg m −3 (standard deviation: 18 kg m −3 ).
Measurements obtained in autumn showed a mean density of 556 kg m −3 (n = 6, standard deviation: 15 kg m −3 ) at the position nearest to the snow surface (about 3.7 m above ground and 1 m below the snow surface), increasing to 578 kg m −3 (n = 6, standard deviation: 8 kg m −3 ) about 2.1 m above ground and to 681 kg m −3 at the lowest sam- The Cryosphere, 12, 385-400, 2018 ple position (n = 1).The resulting bulk density (mean of all three levels) was 606 kg m −3 .This is in accordance with a densification of about 0 to 23 % (9 % for the mean values) in relation to the initial values.

Measured snow volume and mass balance
Snow-height maps calculated from the TLS data of the two snow heaps are presented in Figs.5c, d and 6c, d for Davos and Martell, respectively.The respective average values are shown in Table 4.The maximal height of the Flüela heap was 8.99 m (Martell: 7.61 m) in spring and decreased to 7.86 m (Martell: 6.37 m) in autumn.Mean HS (including sawdust layer) reduced from 4.07 m (Martell: 3.33 m) to 3.15 m (Martell: 2.18 m) and the maximum change in snow height (dHS) was 2.70 m (Martell: 2.32 m).
Red colours in Fig. 5e indicate that dHS was highest at the crown of the snow pile at Flüela.The rather steep crest that was present in spring clearly levelled during the sum- mer, resulting in the formation of a small plateau, a few metres wide, with several small depressions.Contrary to Flüela, the geometry of the Martell heap was much flatter and no such remarkable changes in surface characteristics could be detected, such as the formation of a plateau (Fig. 6a, b).Nevertheless, red colours near the edge of the Martell pile point at higher snow ablation in these locations (Fig. 6e).These higher dHS might be attributed to increased melt caused by the proximity to the paved road (Fig. 2): it can be expected that the road heated much more strongly due to low albedo and thermal properties and thus contributed additional energy by lateral advection of heat towards the snow heap (Mott et al., 2013(Mott et al., , 2015(Mott et al., , 2017) ) or by long-wave radiation.
In contrast, the blue sections visible at the edges of the piles in Figs.5c and 6c represent areas in which height increased.This is, however, an artefact that can be explained by the relocation of sawdust, caused by gravitational and manmade relocation and by the fresh snow at Martell.Note that www.the-cryosphere.net/12/385/2018/The Cryosphere, 12, 385-400, 2018   HS and dHS shown in Fig. 6c, d and e were not corrected for the remaining snow in the depot and for fresh snow covering the surface during the survey in autumn.The reason is that only a few measurements or estimations of the correction values were available that might well be used for reasonable adjustment of mean values or total snow volumes but are inadequate for correction in high spatial resolution.Complementary to overall lower HS at the edges, this fact contributes to the observation of negative HS at the edges of the heap (Fig. 6c, d) and negative dHS in Fig. 6e.
Apart from this no eye-catching patterns such as a possible influence of exposition are visible.This impression is confirmed by Fig. 7a and c, which show box plots relating dHS to aspect.Figure 7b, however, indicates a negative relationship between dHS and the slope at Flüela.Higher ablation seems to happen in flatter places.The flatter (upper) section of the heap is expected to have higher wind speeds and thus more ablation (Fig. 10).Also, the sky-view factor (the amount of sky in the hemisphere above the point of observation versus the amount in the surrounding terrain) is highest there, allowing for a more efficient radiative cooling to the cold sky.For Martell (Fig. 7d), no clear relationship between dHS and the slope can be detected.
The calculated spring snow volume of the Flüela snow pile, including covering material, was 6862 m 3 (Martell: 7138 m 3 ) and decreased to 5307 m 3 (Martell: 4820 m 3 ).This corresponds to a reduction of 23 % (Martell: 32 %).If we relate volume reduction to the initial volume of snow after removing the estimated amount of covering material (Flüela: 830 m 3 , Martell: 860 m 3 ), these values increase to 26 % for Davos and 37 % for Martell.
Nevertheless, it needs to be considered that volume loss and dHS are not only attributed to snow ablation but also to densification of the snow by settling.The density measure- ments described above point to a contribution of about 9 %.This would mean that the effective snow ablation reduces to 17 % (28 %) or -in other words -that 72 to 83 % of the snow mass that had been covered in spring could be conserved over the summer.Note that some snow will additionally be lost during de-covering and distribution of the snow to the ski tracks.
Results obtained from the measurements are compared to results from the SNOWPACK simulations in the next section.It is, however, not meaningful to relate model results calculated for a single point to volume changes for the entire snow pile.We therefore chose to relate model results to the respective coordinates (and settings of the model run), which are the locations of maximal HS.A quadratic area of 1 m 2 surrounding this point was used for the calculation of measured dHS.At these areas the pile shrank by 1.35 m (standard deviation = 2.5 cm) at Flüela and by 1.25 m (standard deviation = 2 cm) at Martell, corresponding to a relative reduction of 15.1 ± 0.3 % (Martell: 16.6 ± 0.2 %).

Model results
Figure 8 illustrates the temporal evolution of HS and density of the 40 cm sawdust-covered snow heap in the Flüela valley.The pronounced yellow line parallel to the surface reflects the lower boundary of the sawdust layer.The initial density of this layer (323 kg m −3 ) rapidly increased in the first 3 weeks due to rising water content caused by rain (initial value 20 %).Afterwards density and water content remained at relatively constant levels of about 580 kg m −3 and 45 % respectively, meaning that the sawdust was always wet.Higher values cannot be reached due to a threshold in the model settings of the layer as described in Sect.2.4.Liquid water fraction of the snow did not change considerably in relation to the initial 3 %.This is attributed to the maximum waterstoring capacity of a snow texture which is implemented in SNOWPACK as derived from Lütschg (2005).
The height of the heap decreased steadily from 9 to 7.60 m.This corresponds to a relative decrease in HS of 15.6 % and is in accordance with the result obtained from the measurements at the point of maximal HS.The results of the reference simulation are very similar for the Martell snow pile (not shown): height decreases by 17.1 % from 7.60 to 6.30 m, which is again an excellent match with the measurements.
Figure 8 also shows densification.It is illustrated how snow is compacted from top to bottom and over the course of time.While simulations were initiated with a constant density of 553 kg m −3 , density at 8 October ranged from 570 kg m −3 at the top to 617 kg m −3 at the bottom of the profile.A comparison with the measurements presented in Sect.2.2.2 shows that the model overestimated snow density in the upper sections of the profile and underestimated it near to the bottom.The simulated bulk density of the entire profile was 600 kg m −3 .Hence, snow densified by 8 % during the summer season.These values agree well with the mean density calculated from the measurements (606 kg m −3 , 9 %).Surface temperature (TS) varied between 0 and 33 • C for the reference scenario and showed a diurnal variation in the range of 5 to 20 • C (Fig. 9).

Sensitivity study of depth of covering layer and meteorology
Figure 10a illustrates the influence of the depth of the sawdust layer.While the entire snow heap melted by mid-September when no covering layer was used, dHS was reduced to 15.5 % (dSWE: 7.2 %) for the 40 cm deep sawdust layer (reference) and finally to 11.1 % (dSWE: 1.5 %) for a 60 cm deep cover.The curves are characterized by an exponential decrease and are clearly flattening out at a thickness of about 40 cm.A reason for the flattening is that thicker layers decrease temperature gradients between the snow surface and atmosphere.Doubling the height of sawdust, for example, reduces the gradient by a factor of 2. As a result, energy input to the snow diminishes as more energy can be stored in the larger sawdust volume and thermal conduction to the snow decelerates.In addition, thicker covers add to dampening effects on TS as illustrated in Fig. 9. Amplitude and extrema, especially minima, are clearly enhanced for shallower covering layers, resulting in higher TS during the day and remarkably cooler TS at night.Heating of the sawdust during the day and cooling during the night appear to be delayed by a few hours.This results in regular changes in the direction of the energy fluxes at the sawdust surface.On average, however, TS (13.8 • C for 40 cm sawdust) clearly appeared to be warmer than TA (11.3 • C), resulting in net negative heat fluxes of sensible heat and long-wave radiation.This mean temperature difference between surface and air clearly decreases for smaller sawdust heights (11.8 • C for 20 cm), consequently diminishing cooling effects by sensible flux and long-wave emission (Table 5).Once a cover layer is thick enough to prevent the surface temperature from dropping significantly below the air temperature, the insulation of a snow heap works well.It can then be assured that it has enough capacity (mass, specific heat) to store the energy during the day while not conducting much to the snow and releasing en-  ergy efficiently at night.Increasing layer thickness exceeding a certain limit of about 40 cm shows only minor improvements in the insulation (Fig. 10a).Temperature gradients in the deeper parts of the cover are then no longer affected by daily changes.Moreover, the relative increase in the layer decreases with thicker cover.In summary, Fig. 10a indicates that a depth of about 30 to 40 cm is required to reach volume savings of around 80 %, while the effect of additional sawdust is minor.The contribution of the different terms of the energy balance at the surface of the heap is shown in Fig. 11a and Table 5, where positive values designate a flux towards the heap (energy source) and negative fluxes a direction to the atmosphere (energy sink).Note that the terms presented here are net values cumulated over the entire simulation period.The huge effect of the covering layer on snow ablation is best illustrated when comparing the energy balance of the reference simulation (40 cm, fifth column in Fig. 11a) to a run without any cover (first column).In total, energy available for ablation is nearly 12 times higher for the simulation without cover.Short-wave radiation is by far the largest source of energy (Table 5, Fig. 11).Due to the much higher albedo of pure snow, net short-wave radiation reduces by about onefifth without sawdust.In contrast, long-wave radiation acts as an energy sink for both simulations but is nearly 13 times higher when a 40 cm sawdust cover is present.All other terms of the surface energy balance, namely net sensible and latent heat fluxes and precipitation, change in sign: while they contribute to melting without cover, they remarkably cool the sawdust-covered snow heap and therefore limit snow ablation.The greatest effect is clearly attributed to long-wave emission, but sensible and latent heat fluxes also contribute substantially.Precipitation only plays a minor role.
The cooling effect of the latent flux is mainly attributed to evaporation at the moist surface of the sawdust.As this layer remained wet for the entire summer and for all simulations shown in Fig. 10a, sawdust thickness appears less relevant to the magnitude of this flux.Finally, Fig. 11a points out that rain, even if only marginally, adds to cooling of the heap, at least when the covering layer exceeded 20 cm.This can again be explained by temperature differences between the warmer sawdust layer and the colder rain.
In addition effects of changed atmospheric forcing, namely WS (Figs. 10b, 11b) TA (Figs. 10c, 11c) and P (Figs. 10d, 11d) have been analysed.WS is an important parameter that especially affects turbulent fluxes (Schlögl et al., 2016).Figure 10b shows that higher WS slightly altered snow ablation.Triplicating WS, for example, resulted in an increase in dHS from 15.5 to 17.6 %.Setting the wind to 0.3 m s −1 reduced dHS to 15.4 %.As expected, rising TA also altered snow ablation (Fig. 10c).For example, an increase in temperature of 5 • C added 25 cm or 3 % to dHS.On the contrary, additional rain did not affect snow ablation significantly (Fig. 10d).Even five times the amount of rain did not change melt considerably (23 % for dHS).The reason is the wetness of the sawdust layer that never dried, as described earlier.Switching P completely, however, increases the snow loss considerably.Smaller difference between HS and SWE for the simulation without P can be explained by less densification of the snow due to reduced moisture content in the snowpack.
In conclusion Figs. 10 and 11 clearly show the effects of the covering layer on the one hand and of the atmospheric forcing on the other hand.This underlines the high impact of sawdust thickness on energy available for snowmelt.Higher WS and warmer TA also altered snowmelt significantly while additional P did not play a decisive role.In general, we rate the accuracy of the snow volumes calculated from TLS data as very high.As described before (Sect.2.2.1), short measurement distances, convenient scan angles, high point densities, and overlapping multiple-scan positions provide favourable conditions for highly accurate measurements.Based on operating experiences with similar settings we estimate the vertical accuracy of the TLS measurements to about 1 cm.Nevertheless, scan shadows still caused some data gaps, especially at the crown of the heaps for the two surveys in autumn (Figs. 5b,6b).These data gaps were caused by a rougher surface such as local depressions.
They therefore had to be closed by linear interpolation, introducing some uncertainty.The extent of the gaps is, however, limited to a few square metres, meaning that the effect on the total mass balance can be rated as marginal (below 0.1 % if we assume an area of 10 m 2 and a mean deviation of 10 cm).Another source of potential error is introduced by the lack of (accurate) bare-ground elevation models.For Flüela no such model could be monitored but the flat and only slightly sloped ground area allowed for a good approximation with a sloped plane defined by the margins of the snow heap.For Martell, a bare-ground elevation model was measured after most of the snow had been distributed in autumn.The remaining snow in the depot, however, could only be estimated based on visual impression and the rating of the local expert.As described in Sect.2.2.1 we assumed snow volume to be in the range of 600 to 1000 m 3 and used 800 m 3 for the corrections.Applying the maximum or minimum estimates would reduce/alter relative snow volume loss by 1 to 2 %.Nonetheless, this correction only affects snow volumes.Snow volume changes or dHS calculations do not require a bare-ground elevation model and are therefore not affected.Finally, thickness of sawdust and fresh snow were obtained from a limited number of probe measurements (see Sect. 2.2.1).Nevertheless, a bias of a few centimetres would be small in relation to the large volume and HS of the snow heaps.
When analysing SWE instead of HS or volume, uncertainty in snow density measurements must also be considered.We showed that the range of snow densities was considerable (541 to 681 kg m −3 ), with higher densities near to the ground.Throughout the storage period a load-and timedependent densification must be assumed due to creep and wet-snow metamorphism.As in spring only densities near to the surface could be measured, an adequate initial density profile with an expected increase with depth could not be captured.This limits the ability to assess the temporal evolu-tion of densification based on measurement but also from the simulation results as the initial snow profile had to be defined with a constant density.

Comparison of the sites
As described earlier, the results of the snow-volume measurements suggested large differences in terms of snow-volume loss between the two sites.While only 22.6 % of the volume disappeared at Flüela, the decrease was more than 10 % larger at Martell (Table 4).Considering the large similarity in terms of initial volume, surface area (Table 4) and also meteorological conditions (Figs. 3,4) this huge difference appears surprising.Possible causes are the different properties of the covering materials (see Sect. 2.1).If we, however, consider the very similar dHS at the position of maximal HS (15 % at Flüela and 16 % at Martell) and the fact that the simulations reproduced this amount very well, this explanation appears insufficient.A second possible reason is a potential warming effect of the black-paved road at Martell (Fig. 2), resulting in a lateral advection of heat as described earlier (Mott et al., 2013(Mott et al., , 2015(Mott et al., , 2017)).The relatively high dHS near to the edges of the heap (Fig. 6) supports this hypothesis.Other micro-meteorological characteristics, such as deviating local wind fields and their implication for energy balance, might be present and could also play a role.Unfortunately, due to a lack of measurement stations directly on the heaps, such effects could not be detected and therefore not proved.Another explanation for increased melt could be insufficient covering of the heap at the edges.Areas with shallow cover were partly visible during our surveys.

SNOWPACK
The excellent agreement of measurements and model results indicates that SNOWPACK is very capable of reproducing dHS of sawdust-covered snow heaps at the point scale.A direct transferability of these results to total mass loss of the heaps in the three-dimensional space requires caution.As shown in Sect.3.3, the total mass loss was much higher than the loss measured and simulated at the point scale.Spatial effects such as variability of radiation and wind and lateral effects at the edges of the heap would need to be considered to simulate total mass losses.Distributed models such as Alpine3D (Lehning et al., 2006) could be used for this task but would also require more distributed information on meteorological forcing, which was not available in this study.
In principal, it needs to be considered that SNOWPACK incorporates some simplifications that might influence simulation results.The sawdust cover was modelled as a homogenous layer.Potential internal heterogeneity such as variable temperature or water content in the layer and its implication to energy balance is not considered.Furthermore, model outcomes are strongly influenced by the accuracy of the input parameters.For example, small-scale variability in meteoro-logical conditions might be present at the test field but not be represented in the recorded meteorological measurements.Moreover, several initial settings such as sawdust properties and the initial state of the snowpack are only based on estimations or limited measurements (see Sect. 2.2).However, sensitivity runs (not shown) that had been performed for the most important parameters (thermal conductivity, water content and texture of sawdust) showed no significant effect on the mass balance.Initial snowpack characteristics (density, water content) were more significant but still only marginally changed the final results.

Conclusions and outlook
A detailed study on snow farming for touristic applications has been presented.Mass balances of two snow heaps covered with a 40 cm thick layer of sawdust and chipped wood respectively have been calculated from repeated TLS surveys in 2015.More than 75 % of the snow volume of a 7000 m 3 snow pile was able to be conserved in the Flüela valley near Davos, Switzerland (Table 4).At the Martell heap, only twothirds of the snow remained in autumn, even though the parameters and conditions were quite similar at both study sites.We assume that this reduced performance is attributed to heat advection from a surrounding paved road or from insufficient cover thickness at some locations on the heap.Moreover, we applied the physically based snow cover model SNOWPACK to simulate snow ablation at the point scale.A comparison of measurements and simulation results showed excellent agreement: at Flüela a dHS of 15.1 % (Martell: 16.6 %) was measured at the point of maximal snow height.Simulations of relative losses (15.6 % for Flüela and 17.1 % for Martell) were only marginally higher.In summary, the magnitude of these results is well in line with operating experiences of different snow-farming sites (12-50 %).
Snowpack simulations were also applied to analyse the contribution of the different terms of the energy balance to snow ablation (Fig. 11).It could be shown that short-wave radiation was by far the most important source of energy.Sensible and latent heat fluxes also contribute to melt if the snow heap is not protected by sawdust.The presence of such a covering layer, however, led to an inversion of the seasonal net fluxes of sensible and latent heat now contributing to cooling of the snow heap.The largest cooling effect was attributed to long-wave emission at the surface of the sawdust.This insulting layer absorbs the solar energy during the day but strongly limits conduction to the snow.At night the absorbed energy is then emitted by long-wave radiation.About two-thirds of the net energy input by solar radiation is compensated for by long-wave cooling.An additional 18 % is balanced by net sensible heat and 6 % by latent heat fluxes.In summary, only 9 % of the net short-wave energy input remained available for snow ablation.The high amount of snow conserved over the summer is therefore attributed to these cooling and insulating www.the-cryosphere.net/12/385/2018/The Cryosphere, 12, 385-400, 2018 effects of the covering layer.Moreover, the high influence of the thickness of the sawdust was evident from the simulations.The larger the covering layer, the smaller the temperature gradient and therefore the smaller the heat flux into the snow.This increased insulation effect finally results in the higher amount of snow that can be conserved.The excellent insulation of sawdust is primarily given by its high density and heat capacity rather than by its heat conductivity (factor of 2 higher than the conductivity of the polypropyleneinsulation plate).This enables the sawdust to absorb large amounts of heat that is re-emitted at night.Still, significant costs of sawdust and additional workload need to be considered when deciding the layer thickness: 40 cm seems a good compromise.
Finally, the effects of varying meteorological conditions, namely WS, TA and P , have been investigated.As expected, it is shown that higher TA and WS enhance snow ablation.Nevertheless, a resulting increase of a few percent appears small in relation to the effect of covering layer thickness.This finding suggests that snow farming might also be feasible under much warmer climatic conditions, indicating that it might also be applicable for lower situated, subalpine ski resorts and communities.Additional precipitation did not play a significant role, as the sawdust remained wet during the entire summer for all scenarios with rain.Switching off rain completely, however, led to a clear increase in snow ablation.As a consequence, irrigation of snow heaps, as suggested by some practitioners, seems unnecessary as long as a certain amount of rain is available.
In conclusion snow farming appears to be an appropriate method for the allocation of a basic snow offer in autumn.However, operating costs and space requirements considerably limit the amount of snow that can be stored.Operational costs have to be evaluated for each snow-farming project, specifically considering the applied technical and logistical solutions.For example, in Davos 15 CHF m −3 snow was estimated for the first snow-farming project in 2008.Until 2016 these costs could be strongly reduced to about 9 CHF m −3 due to larger volumes of stored snow and improved infrastructure and work flow.Investments in structural measures at the storage location are not considered in this calculation.Two-thirds of the expenses were caused by the distribution of the snow along the cross-country track followed by the removal of the sawdust (14 %) and material costs for sawdust (10 %), assuming a 5-year operational lifetime (N.Gruber and W. Putzi, personal communication, 2017).Generally, it can be stated that snow production costs are minor compared to covering costs, especially distribution costs, and therefore cost-efficient snow farming remains limited to applications for which machine-made snow production is impossible or too risky because of weather uncertainties.
SNOWPACK proved to be an appropriate simulation tool and can easily be applied to study the suitability of potential snow-farming sites, to asses the performance of different covering materials (see also Olefs and Lehning, 2010) or to simulate influences on changing climatological conditions on snow farming.The only prerequisite is the availability of adequate meteorological input data.Spatially distributed models such as Alpine3D could be used to investigate spatial effects.The application of multiple different snow models could be interesting when assessing model sensitivity and uncertainty.Future research might also include detailed surveys of other snow-farming projects, possibly at a higher temporal resolution or the investigations of mass losses during different work steps, such as shaping or de-covering snow heaps.Investigations of different types of covering material are currently running in Scandinavia and will contribute to knowledge on best practices of snow farming.

Figure 2 .
Figure 2. Martell snow heap covered with sawdust and wood chips in spring 2015.The Riegl VZ6000 laser scanner can be seen in the foreground.

*
3 • C at Windtunnel and 10.6 • C at Hintermartell (19 May to 27 October).The mean of corrected WS at WindtunnelTable 3. Parameters and settings (offset or factor) applied to the sensitivity studyConstant minimum wind speed of 0.3 m s −1 as hard-coded in SNOWPACK to avoid model instability due to zero surface fluxes of sensible and latent heat.was0.86 m s −1 and the measured WS at Hintermartell (HM) was 1.3 m s −1 .Precipitation cumulated to 543 mm (HM: 534 mm).Net short-wave radiation summed to 1476 MJ m −2 (HM: 1213 MJ m −2 ) and net long-wave radiation was calculated at −992 MJ m −2 (HM: −596 MJ m −2 ) for a heap covered with 40 cm of sawdust.Lower temperatures and irradiation are mainly attributed to the survey period at Martell 3 weeks later.This is confirmed by a comparison of data for an identical period (1 April to 30 September 2015) that demonstrated very similar climatic conditions at the two sites (WT: 10.1 • C, 514 mm; HM: 10.2 • C, 514 mm).Long-term records of temperature and precipitation data from the nearby weather stations DAV and Stausee Zufritt (SZ) were analysed to rate the summer of 2015 in relation to climatic mean values.At DAV the mean temperature of the summer half of the year (1 April to 30 September) 2015 was 10.1 • C (SZ: 9.3 • C) and therefore 0.8 • C (SZ: 0.4 • C) warmer than the mean of the same period of the last 15 years.

Figure 3 .
Figure 3. Meteorological input data aggregated to daily values for the Flüela site.Air temperature (TA) and precipitation are shown in the upper panel, wind speed (WS) and relative humidity (RH) in the middle and net short-wave (SWR) and net long-wave (LWR) radiation in the bottom panel.Measurements were obtained from the stations WT (TA, WS, RH) and DAV (SWR, LWR, P ).

Figure 4 .
Figure 4. Meteorological input data aggregated to daily values for the Martell site.Air temperature (TA) and precipitation are shown in the upper panel, wind speed (WS) and relative humidity (RH) in the middle and net short-wave (SWR) and net long-wave (LWR) radiation in the bottom panel.Measurements were obtained from the station HM, and LWR was parametrized as described above.

Figure 5 .
Figure 5. Flüela snow heap: hillshade images from (a) 29 April and (b) 8 October, snow heights HS (c, d) and snow height change dHS (e).The location of the maximum HS on 29 April is indicated in panels (c) and (d).

Figure 6 .
Figure 6.Martell snow heap: hillshade images from (a) 19 May and (b) 28 October, snow depth HS (c, d) and snow height change dHS (e).Displayed are raw data that were not corrected for the remaining snow in the depot (c, d) and fresh snow at the surface (d, e).The location of the maximum HS on 19 May is indicated in panels (c) and (d).

Figure 7 .
Figure 7. Box plots relating dHS to aspect (a, c) and slope (b, d) of the spring survey for Flüela (a, b) and Martell (c, d).

Figure 8 .
Figure 8. Simulated temporal evolution of snow height and snow density at the Flüela snow heap.

Figure 9 .
Figure 9. Air temperature (red) and surface temperature simulated for the Fluela snow heap with 40 cm (solid blue) and 20 cm (dashed blue) sawdust layer over 4 days in August.

Figure 10 .
Figure 10.Relative loss of snow height (HS) and snow water equivalent (SWE) calculated by SNOWPACK simulations with different thicknesses of the sawdust layer (a), different wind speeds (b), air temperatures (c) and precipitation (d).

Table 5 .
Net energy fluxes at the surface of the Flüela snow heap summed for the entire simulation period without and with a 20 and 40 cm sawdust layer.

Figure 11 .
Figure 11.Most important terms of the energy balance of the Flüela snow heap and their contribution to melt for the different model runs of the sensitivity study: (a) depth of covering layer, (b) wind speed, (c) air temperature, (d) precipitation.
quality of the measurements and resulting error estimates

Table 2 .
Rinderer, 2009)ies of the sawdust cover layer.Note that SNOWPACK characterizes soil and other layers with the volume fractions of the dry material, the void and the water fractions.The properties of the dry material are partly based on measurements and partly on experience (e.g.Rinderer, 2009).

Table 4 .
Geometric properties and snow depth statistics of the Flüela and the Martell snow heap as calculated from the TLS surveys.HS mean is the mean snow height of the heap, HS max the maximum snow height, HS SD the standard deviation of snow height, dHS max the maximum snow height change and dHS at HS max the snow height change at the location of HS max .