Identifying airborne snow metamorphism with stable water isotopes

Wind-blown snow is a frequent phenomenon in high-elevation and polar regions which impacts the surface energy and mass balance of these areas. Loose surface snow gets eroded and transported by wind, which influences the snow particles' physical properties (size, shape, optical properties) that determine the characteristics of the emerging wind-impacted snowpack layer. During airborne snow transport, the governing processes happen on the micro-scale while the particles are transported over long distances. The unfolding processes and the evolution of the particles' physical properties are thus difficult to observe in situ. Here, we used cold-laboratory ring wind tunnel experiments as an interim solution to study the governing processes during airborne snow transport with stable water isotopes as tracers for these micro-scale processes. Repeated analysis of airborne-sampled snow by micro-computed tomography (µCT) documented a growing and rounding of snow particles with transport time, with a concurrent decrease in specific surface area. Stable water isotope analysis of airborne snow and water vapour allowed us to attribute this evolution to the process of airborne snow metamorphism. The changes observed in the snow isotopic composition showed a clear isotopic signature of metamorphic deposition, which requires particle–air temperature gradients. These results question the validity of the thermal-equilibrium assumption between particles and air inside the saltation layer of wind-blown snow events, where the conditions are similar to the ones found in the wind tunnel. Our results thus refine the understanding of the governing processes in the saltation layer and suggest that the snow's isotopic composition can inform on local wind-blown snow events as the original snow isotope signal gets overprinted by airborne snow metamorphism. Within transport times of 3 h, we observed changes in the isotope signal of airborne snow of up to +1.47 ‰ in δ¹⁸O, ±5.7 ‰ in δD, and −6.1 ‰ in d-excess. Thus, airborne snow metamorphism has the potential to influence the climate signal stored in snow and ice core stable water isotope records.

1 Introduction

The phenomena of blowing and drifting snow are frequently observed in snow-covered territory such as polar and alpine regions (Wagner et al., 2022; Walter et al., 2020; Yang et al., 2023). Previously accumulated snow is mobilised and redistributed by wind, which affects both the surface energy and mass balance of the region and feeds back on snow accumulation patterns (Agosta et al., 2019; Lenaerts et al., 2012; Groot Zwaaftink et al., 2011). In the context of intense climate change in the cold climate zones around the globe, it is thus critical to understand the underlying physical processes of wind-blown snow events to accurately model and predict future snow accumulation changes in these regions and their feedback on the climate system.

Wind-driven aeolian transport of snow can be classified into drifting snow close to the ground (∼ 2 m) and blowing snow further aloft. The dominant snow particle transport modes are creep and saltation during drifting-snow events and suspension during blowing-snow events (Bagnold, 1971; Melo et al., 2022). Both transport regimes can coexist, which motivates the terminology of a saltation and suspension layer in aeolian snow clouds (Nemoto and Nishimura, 2004).

The apparent effect of wind-blown snow is the redistribution of snow mass, which, depending on the surrounding topography and wind regime, can form snow surface structures such as ripples, dunes (Amory et al., 2016; Filhol and Sturm, 2015; Sommer et al., 2018a), and cornices (Yu et al., 2023). On glaciers and ice sheets, snow redistribution also directly affects the surface mass balance when snow is blown off the ice margins, although this contribution is small (Palm et al., 2017). More importantly, wind-blown snow contributes negatively to the surface mass balance by enhancing sublimation of snow (Gerber et al., 2023; Palm et al., 2017; Sigmund et al., 2021). Concurrently, enhanced sublimation during wind-blown snow events leads to enhanced energy loss through a near-surface latent heat flux, which affects the near-surface energy balance and can, in turn, negatively feedback or self-inhibit the sublimation process by lowering the air temperature and increasing the relative humidity of the surrounding air (Amory and Kittel, 2019). Apart from these direct impact pathways on the surface energy balance during actual aeolian snow transport, wind-blown snow also indirectly affects the surface energy balance after the airborne snow has settled on the ground. This post-event influence is due to the modification of the particle size and shape of airborne snow particles during wind-blown snow events. We will use the term physical properties of particles (PPPs) in analogy to the physical properties in ice core studies when referring to snow particle shape, size, and morphology. After settling, the PPPs define energy-balance-relevant characteristics of the resulting snowpack such as albedo (Domine et al., 2006), emissivity (Hori et al., 2013), density, and specific surface area (SSA) (Proksch et al., 2015; Schleef et al., 2014b). Previous studies (Sommer et al., 2017; Walter et al., 2024) have shown that wind-driven aeolian transport of snow significantly affects the microstructure of the ultimately deposited surface snow; i.e. density increases and SSA decreases are intensified under the influence of wind. A resulting wind-blown snow layer is typically characterised by small, broken, or abraded; closely packed; and well-sintered particles (Fierz et al., 2009).

Mechanisms that can control the modification of the PPPs and thus the resulting snowpack characteristics can be grouped into mechanical processes such as particle fragmentation, abrasion (Clifton et al., 2006; Comola et al., 2017; Gromke et al., 2014), and aggregation (Lo and Passarelli, 1982) or into thermodynamic (hereafter metamorphic) processes such as sublimation at snow particle surfaces (Schmidt, 1982; Thorpe and Mason, 1966) and vapour deposition (re-sublimation) on the suspended snow particle (Sharma et al., 2018; Sigmund et al., 2024; Yamaguchi et al., 2019). In this context, Walter et al. (2024) first introduced the term “airborne snow metamorphism”, which summarises the multiple cycles of sublimation and vapour deposition on the suspended snow particle resulting in modifications of the snow particle size and shape during aeolian particle transport in analogy to metamorphism inside a stationary snowpack (Pinzer et al., 2012; Schleef et al., 2014b). The term snow metamorphism describes the recrystallisation of snow grains in a snowpack that is driven by vapour pressure gradients (Colbeck, 1982). Snow metamorphism typically results in a decrease in SSA (growing of snow grains) and can be associated with density changes in a snowpack (Jafari et al., 2020; Kaempfer and Schneebeli, 2007). Based on the temperature regime of the snowpack, a distinction can be made between isothermal and temperature gradient metamorphism. The dominant physical processes creating the vapour pressure gradients are distinct: under isothermal conditions, the micro-scale curvature effects (Kelvin equation) drive metamorphism (Colbeck, 1980), which are outweighed by macro-scale temperature gradient effects (Clausius–Clapeyron equation) in snowpacks with a temperature gradient (Marbouty, 1980). Temperature gradient metamorphism typically results in higher recrystallisation rates and thus faster SSA decay in snowpacks (Taillandier et al., 2007).

Thus, airborne snow metamorphism was proposed as a driving factor for PPP changes during aeolian transport of snow, yet the relative importance of the different processes involved and their combined effect on the snow microstructure are still unknown due to missing observations. Snowpack or climate models that include wind-blown snow modules thus rely on field-based empirical parameterisations when incorporating wind-blown snow effects on snowpack characteristics because the physical process understanding is limited (Walter et al., 2024).

In situ observations of changes in the PPP during wind-blown snow events are difficult due to the micro-scale nature of these processes and the short timescales on which they evolve while the particles cover large distances on their airborne trajectories. Thus, in situ observations are limited to point measurements within the snow cloud (Nishimura et al., 2014; Nishimura and Nemoto, 2005) or pre- and post-event monitoring of the snowpack characteristics (Walter et al., 2024). Field measurements are the foundation for wind-blown snow parameterisations that are currently implemented in snowpack and climate models (Lehning et al., 2002; Liston et al., 2007; Vionnet et al., 2012). As an approach to circumvent the difficulties of multiple unconstrained environmental variables during in situ field observations, controlled laboratory experiments using straight wind tunnels (Nishimura et al., 2014; Wever et al., 2009) or ring-shaped wind tunnels with infinite fetch (Sommer et al., 2018b; Walter et al., 2024) have been used successfully to systematically validate model parameterisations. Additionally, shadowgraphy imaging based on high-speed camera recordings (Gromke et al., 2014) has been applied in wind tunnel experiments to investigate the vertical distribution of snow particle size and number. However, this method does not allow access to information on the evolution of PPPs and thus on the differentiation between mechanical (fragmentation, aggregation, and abrasion) and metamorphic (sublimation and vapour deposition) drivers. Consequently, means to constrain the underlying governing equations of the PPP evolution during airborne snow metamorphism are still missing due to the nature of the processes unfolding on the difficult-to-observe micro-scale.

As direct observations of the individual airborne snow particles are currently impossible, indirect methods such as the monitoring of proxy variables are key to being able to disclose metamorphic processes happening on the micro-scale in macro-scale variables. For snow processes, such a proxy can be the abundance of stable water isotopologues (hereafter simplified to “water isotopes”), namely H₂¹⁶O, H₂¹⁸O, and ²H¹H¹⁶O, in a sample. Due to their predictable partitioning, i.e. fractionation into the vapour, liquid, and ice phases during phase changes (i.e. water vapour deposition and sublimation or freeze and melt), stable water isotopes are powerful tracers (Ala-aho et al., 2021; Beria et al., 2018; Galewsky et al., 2016) which allow us to identify and quantify the impact of metamorphic processes on the PPPs (Ebner et al., 2017; Harris Stuart et al., 2023). Ideally, stable water isotope measurements comprise all involved phases (in dry-snow regions this, means solid and gas phases) to fully constrain the isotopic fractionation during phase change processes. Stable water isotopes are measured as the abundance ratio R of the rare (i.e. ¹H₂¹⁸O and ²H¹H¹⁶O) to the most abundant (¹H₂¹⁶O) isotopes and are commonly reported as δ values referenced against an international standard (equal to the isotope ratio of the mean ocean water) and are given in per mil (‰) (Craig, 1961b; Mook, 2000):

$$\begin{array}{}\text{(1)}& {\mathit{\delta }}^{\mathrm{18}}\mathrm{O}={\displaystyle \frac{R{{}^{\mathrm{18}}\mathrm{O}}_{\mathrm{sample}}}{R{{}^{\mathrm{18}}\mathrm{O}}_{\mathrm{standard}}}}-\mathrm{1}(\cdot {\mathrm{10}}^{\mathrm{3}}\mathrm{‰}).\end{array}$$

Isotopic fractionation is temperature (T) dependent and empirically well established for equilibrium fractionation processes between vapour and ice (Majoube, 1970; Merlivat and Nief, 1967). For non-equilibrium processes, i.e. processes including a humidity gradient (RH ≠1), the net fractionation between two phases deviates from the pure equilibrium fractionation factor α_eq by a non-equilibrium (sometimes also termed kinetic) fractionation factor α_k that is due to different molecular diffusivities for the different water isotopes (Merlivat, 1978). The effective fractionation factor of a phase change is thus defined as follows:

$$\begin{array}{}\text{(2)}& {\displaystyle \frac{R{{}^{\mathrm{18}}\mathrm{O}}_{\mathrm{ice}}}{R{{}^{\mathrm{18}}\mathrm{O}}_{\mathrm{vapour}}}}={\mathit{\alpha }}_{\mathrm{net}}\left(T,{\mathrm{RH}}_{\mathrm{s}}\right)={\mathit{\alpha }}_{\mathrm{eq}}\left(T\right)\times {\mathit{\alpha }}_{\mathrm{k}}\left({\mathrm{RH}}_{\mathrm{s}}\right),\end{array}$$

where RH_s is the relative humidity normalised to surface saturation and defines the humidity gradient. Usually, α_net>1, which shows that the resulting ice phase is enriched in the rare water isotope compared to the relatively depleted vapour phase; hence, ${{\mathit{\delta }}^{\mathrm{18}}\mathrm{O}}_{\mathrm{ice}}\mathit{>}{{\mathit{\delta }}^{\mathrm{18}}\mathrm{O}}_{\mathrm{vapour}}$. However, deposition under supersaturation conditions (RH_s>1) can lead to strong non-equilibrium fractionation, such that ${{\mathit{\delta }}^{\mathrm{18}}\mathrm{O}}_{\mathrm{ice}}\mathit{<}{{\mathit{\delta }}^{\mathrm{18}}\mathrm{O}}_{\mathrm{vapour}}$ (Jouzel and Merlivat, 1984). Equations (1) and (2) can be written analogously for ¹H²H¹⁶O, with the notation δD being used onwards. From both isotopes, the second-order variable d-excess can be derived as $\text{d-excess}=\mathit{\delta }\mathrm{D}-\mathrm{8}\times {\mathit{\delta }}^{\mathrm{18}}\mathrm{O}$ (Craig, 1961a). The slope of 8 originates from the 8-fold stronger fractionation of ²H¹H¹⁶O compared to H₂¹⁸O under equilibrium conditions. Per definition, this variable will not change under pure equilibrium conditions, and changes in d-excess thus indicate non-equilibrium fractionation processes. Note that, due to the non-linearity of the δ scale, under very depleted conditions, this assumption does not hold (Dütsch et al., 2017). An additional concept to understand and predict isotope behaviour is the analysis of the out-of-equilibrium state, i.e. the disequilibrium in isotopic composition between two phases of water, similarly to the analyses by Aemisegger et al. (2015) or Wahl et al. (2022). For this, the temperature-dependent equilibrium state between two phases is compared against the actual isotopic composition of the respective phase. The difference between the equilibrium and the actual isotopic composition is called disequilibrium, and it describes whether a water entity is more enriched or depleted compared to the equilibrium state. This analysis can inform us about the expected evolution during the mixing of water entities of the original and newly produced phases and can help to reconcile isotope signal behaviour in transient conditions.

This paper aims to characterise airborne snow metamorphism and its effect on the temporal evolution of PPPs during wind-blown snow events for the first time through the lens of stable water isotopes. To this end, we set up cold-laboratory experiments using a ring wind tunnel to simulate wind-blown snow events under controlled conditions and for long transport durations (Sect. 2). We monitored the co-evolution of PPPs and the isotopic composition of both aeolian snow and water vapour inside the ring wind tunnel. From these observations, we deduce the underlying physical processes that are responsible for changes in PPPs and relate them to the ambient atmospheric conditions (Sect. 3). We discuss the significance of our findings for resulting snowpack characteristics and for the modelling and calculation of the energy and mass balance of wind-blown snow events, and we reflect on the implications of our findings for the climate signal imprinted in the snow isotopic composition (Sect. 4).

2 Materials and methods

All experiments were conducted at the cold-laboratory facilities of the WSL Institute for Snow and Avalanche Research SLF in Davos, Switzerland. The experiments were performed in two periods, one from 9 to 13 January 2023 and the other from 11 to 31 May 2023. In the following, an overview of the 19 experiments and the sampling protocol is given, and the experimental setup is described.

2.1 Cold-laboratory wind tunnel experiments

The experiments aimed to simulate wind-blown snow transport in the ring wind tunnel for long transport times while monitoring changes in snow properties and environmental variables to identify the corresponding governing mechanisms. For this, 600 g of fresh snow with known initial PPPs and isotopic composition were introduced to the wind tunnel at an initial wind speed of ∼ 8.5 m s⁻¹, and airborne snow was then repeatedly sampled out of the airstream until no more snow was in suspension.

To keep the experiments' initial conditions as repeatable and reproducible as possible, we used artificially produced snow from SLF's snowmaker (Schleef et al., 2014a), which produces nature-identical snow from tap water. The snowmaker was situated in the same cold-laboratory, but the snowmaker and the wind tunnel were not operated simultaneously to limit heat and moisture sources during the experiment hours. The snow was produced at an air temperature of −20 °C and a water bath temperature of +30 °C, yielding dendritic and hollow column crystals with very low density (45–80 kg m⁻³) and high SSA (50–70 mm⁻¹) (Walter et al., 2024). The produced snow was stored at −20 °C for a maximum of 4 d before it was used for the experiments. Care was taken to equilibrate the snow temperature to the target wind tunnel air temperature (Table 1) for 30–60 min before the start of each experiment. The snow was well mixed, its temperature was measured, and three replicate samples of 10 g were taken before the introduction of the snow to the wind tunnel to estimate the initial isotopic composition. The mean of these three samples serves as the initial-condition snow isotope value for each experiment. The three initial-snow samples agreed to within a range of standard deviations (SDs) of 0.01 ‰–0.66 ‰ for δ¹⁸O, 0.0 ‰–2.8 ‰ for δD, and 0.1 ‰–2.5 ‰ for d-excess, with a mean SD of 0.12 ‰ for δ¹⁸O, 0.7 ‰ for δD, and 0.6 ‰ for d-excess within the set of 19 experiments.

At the start of the experiments, the wind tunnel was closed before the snow was prepared for introduction. After the conditions inside the wind tunnel had reached an equilibrium, the wind generation was initiated quickly, reaching wind speeds >8 m s⁻¹, which is well beyond the threshold wind speed for snow movement and ensures that snow is transported in suspension (He and Ohara, 2017; Yu et al., 2023). Shortly after reaching the maximum wind speed, the snow was added to the wind tunnel in small batches through a port in the ceiling within a short time window of 2–5 min. After the snow introduction to the wind tunnel, a part of the snow was not kept in suspension due to the geometry of the flow field in the wind tunnel (Sommer et al., 2017). Some snow accumulated in the curved section underneath the propeller location. In addition to this early and permanent deposition of snow, a saltation layer with a height of approximately 5–10 cm (Walter et al., 2024; Yu et al., 2023) formed in the wind tunnel, where snow was continuously deposited and re-eroded. Importantly, some part of the snow was permanently in suspension up to heights of 10–40 cm. Throughout the whole experiment and until all airborne snow particles had disappeared (due to sublimation or accumulation), the environmental conditions (temperature (T), wind speed (WS), relative humidity (RH)) and the vapour isotopic composition (δ_v) inside the wind tunnel were monitored using the sensors presented in the next section. Additionally, a small amount of snow (3–10 g) was sampled from the air several times during a single experiment using a Whirl-Pak^® that was tied to a hand net frame. This sampling tool was inserted into the wind tunnel through a hole in the wind tunnel ceiling and was moved up and down in the air stream, covering the whole height of the wind tunnel until enough snow had been captured in the Whirl-Pak^® for later water isotope analysis in the laboratory. During six experiments, part of the so-collected snow was transferred to specific sample holders for the analysis of particle shape, size, and distribution using micro-computed tomography (µCT). Thereafter, the Whirl-Pak^® was closed and stored at −20 °C until further preparation for shipment to the isotope laboratory facility. Once no snow was left in the air column, the experiment was declared to be finished, and the propeller was turned off. Two samples of the accumulated snow were taken right after the wind tunnel was opened: one from the surface layer and one from inside the accumulated snow. Before the start of a new experiment, the sealed wind tunnel was opened to the laboratory air to allow sublimation of vapour depositions on the walls and ceiling, and all the leftover snow was cleaned out.

The temperature regime at which the experiments were conducted was regulated by adjusting the room temperature of the laboratory. During the experiment runs, the air temperature inside the wind tunnel was observed to increase, whereas the room temperature of the cold-laboratory was kept stable (±1 °C) by the laboratory climate control system. The source of this observed warming of the air within the wind tunnel during the experiment can be primarily attributed to the heat generation of the motor and propeller setup. A control run with a sealed but empty wind tunnel was performed to characterise the experiment setup behaviour without added snow. The measurements of this control run are presented in Fig. S1 of the Supplement (SM). It was identified that the precise position of the wind tunnel inside the laboratory with regard to the internal air-conditioning system played a major role in limiting the heating of the air inside the wind tunnel during the experiment. The heating was minimised when the wind tunnel was placed directly opposite the air conditioning system, which allowed for the most effective air circulation in the cold laboratory.

In total, 19 experiments were conducted at different temperature regimes, spanning a temperature range of −20 to −3 °C at the start of the experiments. In two special cases, the laboratory temperature was increased or decreased during the experiment duration, leading to pronounced warming (experiment no. 18, $\mathrm{\Delta }T:+\mathrm{13.5}$ °C) or cooling (experiment no. 12, ΔT: −6 °C) in the wind tunnel. Details on the specific characteristics of each experiment are summarised in Table 1.

Table 1 The characteristics of the 19 ring wind tunnel experiments. During the three experiments in January (^∗), no water vapour isotopes were measured, and the humidity increase (ΔH₂O) is estimated based on the RH sensors (P_w refers to water vapour pressure, and q refers to specific humidity). ΔT lists the change in air temperature inside the wind tunnel during the respective experiment (disc. indicates discarded).

2.2 Instrumentation

2.2.1 Ring wind tunnel and setup

We used an existing ring wind tunnel (Fig. 1; Yu et al., 2023; Sommer et al., 2018b) that was placed in a cold laboratory of the SLF. An outer and inner ring made of wood and plexiglass are placed on wooden planks a few centimetres above the laboratory floor, forming a tunnel with a rectangular cross-section and dimensions of 0.5 m height and ∼ 0.2 m width. The resulting oval ring has outer dimensions of 2.20 m length and 1.20 m width, and the removable lids can be attached to the wooden frame with metal locks. The wind tunnel has an air volume of ∼ 0.5 m³ and was sealed against snow loss with the help of insulation material. However, the tunnel is not fully air-tight, and minor exchange with the laboratory air is inevitable. Wind is generated by an electrically motored two-blade propeller that accelerates the air at a height of approximately 40 cm in the wind tunnel. The wind tunnel is equipped with basic meteorological instrumentation in the observation section of the tunnel, where the lid and walls are made from plexiglass (Fig. 1). This is where WS is recorded by a MiniAir60 (Schildknecht) wind sensor at 37 cm height. HygroClips (HC2-S) (Rotronics) measuring RH and T are installed at three heights (15, 30, 42.5 cm), and an additional Pt100 (K. Mösch AG) T sensor is installed at 3 cm. The T and RH sensors were calibrated in dedicated calibration runs. During these runs, they were placed close together and referenced against their mean readings in a range of −21.5 to −5 °C and 46 % to 91 % to be comparable. For water vapour sampling, we installed a heated copper tube (OD $\mathrm{1}/\mathrm{4}$”, ID 0.19”) that entered the wind tunnel through a hole in the lid at the same observation section. The copper tube was bent such that its opening was oriented with the wind direction, looking downwind at a height of 47 cm in the wind tunnel. The tube then was fed through a connection hole in the laboratory wall and connected to a KNF membrane vacuum pump that was installed in the temperate surroundings outside the cold laboratory. The tube was heated to about 5 °C above the cold-laboratory temperature and insulated with Armaflex insulation material to prevent condensation inside the tube. The flow rate inside the copper tube was set to 200 mL min⁻¹ and controlled by a Sensirion mass flow controller of type SFC5500-10SML. The wind tunnel has four round holes with a 10 cm diameter in each curve; these are closed with metal lids but allow access to the inside of the wind tunnel. Samples of airborne snow were taken through one of these holes.

Figure 1 Wind tunnel setup. In (b), the ring wind tunnel is shown as it was installed in the cold laboratory. A close-up of the observation section is shown in (a), displaying the positions of the sensors that measured the environmental conditions inside the wind tunnel during the experiments. The vapour measurement setup including the cavity ring-down laser spectrometer (CRDS) was installed outside the cold laboratory and is shown schematically in (c). The KNF pump pumped air at a low flow rate of 200 mL min⁻¹ out of the wind tunnel (blue arrows) to the CRDS for continuous real-time analysis of the water vapour and its isotopic composition. During CRDS calibrations, the standard delivery module (SDM) and vaporiser generated a vapour stream of constant humidity and isotopic composition to be measured by the CRDS.

2.2.2 Stable water isotope measurements

During the experiments, we measured the water vapour inside the wind tunnel continuously and sampled drifting and blowing snow periodically for stable water isotope analysis of both water phases.

Isotopic composition of the water vapour

For analysing the isotopic composition of the water vapour, a cavity ring-down laser spectrometer (CRDS) of type Picarro L2130-i was connected to the copper tube outside the laboratory at stable room temperature (∼ 18 °C) and drew in air at the instrument's flow rate of 30 mL min⁻¹ through the connected vaporiser. The total length of the copper tube between the inlet in the wind tunnel and the entry to the vaporiser was 4.6 m for the first five experiments and 6.1 m after repositioning the wind tunnel in the laboratory. The time delay between a signal occurring in the wind tunnel before it was recorded by the CRDS was experimentally determined to be 38 and 50 s, respectively. As we expect no abrupt jumps in the vapour isotopic composition from the processes in the wind tunnel but rather smooth transitions, this lag time is acceptable. The meteorological and water vapour readings were synchronised while accounting for this time lag during data post-processing.

The recorded vapour measurements were calibrated following standard procedures.

The humidity readings of the CRDS were calibrated with a LI-610 Dew Point Generator in the temperature range between 0 and 20 °C, with an emphasis on the range between 0 and 2.5 °C. In total, 14 calibration pulses of constant (maximum SD ±16 ppmv) humidity between 20–30 min were generated with the LI-610 and were calibrated against the respective saturation vapour pressure following Murphy and Koop (2005). As independent calibration validation, all CRDS humidity measurements were compared to the humidity measurements of the highest HygroClip when both sensors were measuring at the same time. Both measurement techniques yielded specific humidity (g kg⁻¹) readings that agree reasonably well, with a slope of 0.98, an offset of 0.04 g kg⁻¹, and an R² of 0.99 at p<0.001.

The water vapour isotope measurements were first corrected for the humidity dependency of the CRDS instrument (Aemisegger et al., 2012) and afterwards were calibrated against the VSMOW-SLAP scale. The instrument-specific humidity dependency curve was established using a standard delivery module (SDM) of Picarro by injecting vapour of known isotopic composition at different humidity levels into the CRDS through the vaporiser of model V1102-i in the week prior to the experiments. Since the range of the expected isotopic composition of water vapour was narrow, only one humidity dependency curve was established for a standard water with isotope values close to the expected values (PCD: δ¹⁸O = −20.87 ‰, δD = −169.5 ‰). This water standard was measured in a humidity range between 500–18000 ppmv, with 16 821 ± 44 ppmv as the reference normalisation humidity level calculated from 10 plateaus. The data points for the humidity dependency calibration curves were calculated from 15–25 min of continuous CRDS data at constant (max ±190 ppmv, mean SD of 91 ppmv) humidity levels, and the calibration curves for both isotopes were established as orthogonal least-square fit of form $y=\frac{a}{x}+bx+c$ to a total of 115 data points (SM Fig. S2) following Weng et al. (2020).

For the normalisation to the VSMOW-SLAP scale, two standard waters (MP: δ¹⁸O = −11.12 ‰, δD = −79.8 ‰; GRIP: δ¹⁸O = −34.14 ‰, δD = −265.0 ‰) bracketing the expected vapour isotopic composition were vaporised and injected into the CRDS by the SDM once a day, first at a high humidity level of 26 000 ppmv for 15 min to reduce memory effects and afterwards at a constant (max SD of 190 ppmv, mean SD of 99 ppmv) humidity level of 16 770 ± 354 ppmv (mean and SD of all 36 accepted calibration pulses). Continuous measurements of 16–25 min were averaged to obtain the calibration point for each day and each standard water. Since no instrument drift was observed during the laboratory campaign, all water vapour measurements were calibrated with the same calibration data set. This is consistent with other campaigns during which negligibly small drifts of the CRDS instrument were observed (Bailey et al., 2023). The calibration data set was generated by calculating a linear least-square fit to all 36 accepted calibration pulses, of which 21 were of standard water MP and 15 of standard water GRIP, obtained from calibration runs on 13 d. The water standards for the SDM had to be refilled once during the campaign, which resulted in small changes in the standard waters' isotopic compositions that were accounted for during the normalisation. All used sets of standard water were analysed for their isotopic composition at the hydrology laboratory of the University of Freiburg. We estimate the uncertainty of the standard water from replicate samples taken a minimum of 5 d apart from the SDM, which had an absolute offset of ≤0.03 ‰ for δ¹⁸O and ≤0.3 ‰ for δD. The humidity-corrected and calibrated 1 s vapour measurements were resampled to 3 min averages. The uncertainty associated with the precision of the vapour isotope data is estimated as the standard error of the mean of the 3 min averages.

The snow samples were melted in the air-tight closed Whirl-Pak^®, filtered with a 0.2 µm filter (VWR, nylon membrane, 25 mm), and 1.5 mL of the sample was transferred to 2 mL glass vials and sealed. The samples were sent to and analysed at the stable water isotope laboratory of the Paul Scherrer Institute in Villigen, Switzerland, following established laboratory procedures with an uncertainty of the individual snow sample isotopic composition of <0.1 ‰ for δ¹⁸O and <0.5 ‰ for δD (Avak et al., 2019), which yields an uncertainty of 0.9 ‰ on d-excess.

2.2.3 Micro-CT scanning

Micro-computed tomography (µCT) 3D microstructure measurements of the airborne snow samples were performed to obtain information on the evolution of the SSA and the particle size distributions of the snow during the experiments. A small portion (≈ 1 cm²) of the periodically sampled snow was filled in a cylindrical µCT sample holder with a diameter of 15 mm. The measurements were performed using a Scanco^® CT-40 scanner (e.g. Heggli et al., 2011; Pinzer and Schneebeli, 2009). A small sample holder was chosen to obtain a high 3D microstructure resolution (voxel size of 8 µm) and thus derive accurate measurements of the SSA and particle size distributions. For the binary segmentation, the energy-based segmentation procedure presented by Hagenmuller et al. (2013) was used. A marching-cubes approach (Hagenmuller et al., 2016) was used for calculating the SSA from the segmented µCT images. A metric for the particle size distribution can be derived from the µCT images by filling the 3D ice matrix with inscribed spheres of different diameters d (Hildebrand and Rüegsegger, 1997).

Table 2 Statistics of the observed changes (Δ) in the isotopic composition of snow (δ_snow) and water vapour (δ_vapour) during the wind tunnel experiments. Only significant changes are used for the calculations of range, mean, and ±1 standard deviation (SD). For Δδ_snow, changes in the first airborne and the maximum observed isotope changes are given (1st $\mathrm{\Delta }{\mathit{\delta }}_{\mathrm{snow}}^{\ast }={\mathit{\delta }}_{\mathrm{snow},\mathrm{1}\mathrm{st}}^{\ast }-\stackrel{\mathrm{‾}}{{\mathit{\delta }}_{\mathrm{snow},\mathrm{initial}}^{\ast }}$; max $\mathrm{\Delta }{\mathit{\delta }}_{\mathrm{snow}}^{\ast }={\mathit{\delta }}_{\mathrm{snow},j}^{\ast }-\stackrel{\mathrm{‾}}{{\mathit{\delta }}_{\mathrm{snow},\mathrm{initial}}^{\ast }}$, where j refers to the airborne sample with largest deviation from initial snow). For Δδ_vapour, the number of experiments with high and low isotope signal changes during and following snow introduction are listed and described further in Sect. 3.2.2.

3 Results

In the following, we present the evolution of the PPPs during the wind-blown snow experiments (Sect. 3.1), followed by the results of the evolution of the snow and vapour isotope signals after calibration and post-processing (Sect. 3.2).

3.1 The evolution of the snow's physical properties during aeolian transport

During 6 of the 16 experiments conducted in May 2023, the airborne snow samples were analysed for their PPPs through measurements of specific surface area (SSA) and sphere size distribution analysis using a µCT scanner. The results of experiment no. 10 had to be discarded as the snow was stored at −2 °C for too long before the samples were measured, which led to distorted µCT results. In the five remaining experiments, the SSA of the last snow sample was reduced to 35 %–70 % of the initial SSA after being transported in the wind tunnel for 1.25–2.5 h. The observed SSA decrease follows an exponential decay as a function of transport time (Eq. 3, Fig. 2a) in agreement with previous SSA decay observations during metamorphism (Cabanes et al., 2003):

$$\begin{array}{}\text{(3)}& \mathrm{SSA}\left(t\right)=({\mathrm{SSA}}_{\mathrm{0}}-c)e^{-\mathit{\alpha }t}+c,\end{array}$$

with t being the transport duration, α being the decay rate, and c being a parameter for the non-zero asymptote. Neither α nor c show a simple dependency on the temperature regime of the experiment.

The sphere sizes detected in the initial and airborne snow samples can be described with a two-parameter gamma distribution, similarly to the expected particle size distributions in wind-blown snow (Budd, 1966; Clifton et al., 2006; Schmidt, 1982). The evolution of the sphere size distribution in all experiments followed a general coarsening pattern. As shown in Fig. 2b, the initial-snow samples (fresh snowmaker snow) show a distinct peak at ∼ 30–40 µm sphere size, with the largest detected sphere sizes reaching 100–150 µm. With time, both the distribution peaks and maximum values gradually shift towards larger sphere sizes. The smallest sphere sizes detected in the fresh snow (∼ 10 µm) had disappeared in the first airborne samples. However, small sphere samples (∼ 20 µm) are universally detected in all samples taken; i.e. there is only a marginal shift in minimum sphere sizes detected with transport time. Hence, airborne transport leads to a broadening of the distributions. After >80 min of airborne transport, the biggest particles have sphere sizes of >175 µm, reaching 200 µm or larger for the last airborne samples taken in each experiment. By the end of each experiment, the most abundant sphere sizes had shifted to 60–75 µm. Interestingly, the samples taken from the accumulated snow show a difference in the distribution of samples taken from the inner part vs. from the surface layer. The sample from the surface layer (“acc_out”, Fig. 2) shows a similar sphere size distribution compared to the airborne samples towards the end of the experiment, whereas the snow from within the accumulation layer (“acc_in”, Fig. 2) typically shows a distribution more closely resembling the early airborne snow samples. This supports our assumption that snow that is deposited early during the experiment remains in the accumulated snow layer and that snow that is sampled from the airstream can be regarded as “active snow”; i.e. it was transported either in saltation or suspension throughout the entire experiment duration. Note that density measurements of snow collected from an airstream are meaningless and are thus not shown.

Figure 2 The results of the µCT scans of five experiments. Panel (a) displays the decrease in SSA values with transport time in all five experiments at different temperature regimes. Dashed lines are exponential fits of Eq. (3) to the data points. Panels (b)–(f) each show the corresponding sphere size distributions for one experiment that serves as a measure for representative snow particle size. Each distribution is calculated from one snow sample, and the colour corresponds to the time after snow introduction. Note that “acc” samples were scraped off the snow accumulation in the wind tunnel. Samples labelled “acc-out” were taken from the surface layer, and “acc_in” samples were taken from the innermost part of the accumulated snow.

3.2 The isotopic signature of wind-blown snow experiments

To identify the underlying process responsible for the changes in the physical properties during aeolian transport, the evolution in the water isotope signal of the snow was obtained from discrete snow sampling at regular time intervals after the snow introduction to the wind tunnel, and the vapour isotopic composition in the wind tunnel was measured continuously.

3.2.1 Changes in the snow isotopic composition

The snow isotopic composition changed with time during aeolian transport. Figure 3 shows the change in isotopic composition between the initial snow (mean of three samples) and the first airborne sample in panels (a), (b), and (c) and the maximum observed change between initial and airborne snow in panels (d), (e), and (f) for δ¹⁸O, δD, and d-excess, respectively. The statistics of the observed significant changes in snow isotopic composition are given in Table 2.

Just after the snow introduction, the majority of the experiments show a significant enrichment in the snow δ¹⁸O value ($\mathrm{15}/\mathrm{19}$) and a decrease in the d-excess value ($\mathrm{12}/\mathrm{19}$). For δD, only $\mathrm{5}/\mathrm{19}$ experiments show a significant enrichment. In all other experiments, the observed changes in the first airborne snow sample are within the uncertainty limits of the calculated difference and thus are not significant. The snow isotopic composition evolved further during aeolian transport.

The general pattern of maximum observed changes between initial and airborne samples that are significant is a δ¹⁸O enrichment ($\mathrm{15}/\mathrm{19}$) and a d-excess decrease ($\mathrm{16}/\mathrm{19}$). Three experiments show a reverse behaviour of significant δ¹⁸O depletion, and one experiment shows a d-excess increase. For δD, the maximum isotope change recorded is almost evenly distributed between enrichment ($\mathrm{8}/\mathrm{19}$) and depletion ($\mathrm{6}/\mathrm{19}$) cases, and five recorded changes fall within the range of measurement uncertainty. Both the cooling (no. 12, ΔT: −6.3 °C) and the warming (no. 18, ΔT: +13.5 °C) experiment follow the general pattern of δ¹⁸O enrichment and d-excess decrease in the first airborne snow sample, which intensifies during the experiments. For δD, both experiments show a significant enrichment only after the first snow-sampling time. During the warming experiment, 14 airborne snow samples were taken, which allows a robust trend analysis for this experiment. Including the initial mean isotopic composition and accounting for the uncertainties, the trends for snow isotopic evolution have positive slopes (m) for δ¹⁸O (m=0.006 ‰ min⁻¹, p=0.000) and δD (m=0.011 ‰ min⁻¹, p=0.091) and a negative slope for d-excess ($m=-\mathrm{0.049}$ ‰ min⁻¹, p=0.002).

No general, simple dependency between observed isotopic changes and average temperature, net sublimation, transport duration, or initial isotopic composition can be identified. Note that the net increase in humidity correlates positively with the average experiment temperature, as expected by the Clausius–Clapeyron relationship (Fig. 3).

Figure 3 Overview of recorded changes in snow isotopic composition during the 19 different experiments. The upper row shows the change in snow isotopic composition (${\mathit{\delta }}_{\mathrm{snow}}^{\ast }$) between the average initial snow and the first airborne snow sample (1st $\mathrm{\Delta }{\mathit{\delta }}_{\mathrm{snow}}^{\ast }={\mathit{\delta }}_{\mathrm{snow},\mathrm{1}\mathrm{st}}^{\ast }-\stackrel{\mathrm{‾}}{{\mathit{\delta }}_{\mathrm{snow},\mathrm{initial}}^{\ast }}$). The lower row shows the maximum observed change between average initial and airborne snow (max $\mathrm{\Delta }{\mathit{\delta }}_{\mathrm{snow}}^{\ast }={\mathit{\delta }}_{\mathrm{snow},j}^{\ast }-\stackrel{\mathrm{‾}}{{\mathit{\delta }}_{\mathrm{snow},\mathrm{initial}}^{\ast }}$, where j is the airborne sample with the largest deviation from initial snow). Changes in δ¹⁸O are displayed in (a) and (d), changes in δD are shown in (b) and (e), and changes in d-excess are shown in (c) and (f), plotted against the average temperature during the experiments. The colour code in (a), (b), and (c) shows the average initial-snow isotopic composition. The colour code in (d), (e), and (f) shows the change in specific humidity inside the wind tunnel due to the sublimation of snow until the sampling of the corresponding snow sample. The warming (red cross) and cooling (black cross) experiments are excluded from the colour coding because of their extreme values, which are given in Table 1. Circular grey data points are the experiments conducted in January, in which the specific humidity increase was calculated from the HygroClip sensors in the absence of the CRDS measurements. The marker size indicates the total transport duration in the range of 49 to 180 min, with bigger markers representing longer transport times. Error bars account for the total error calculated from the standard error of the initial mean and the measurement uncertainty of the individual airborne snow sample.

The snow samples taken from the accumulated snow (not shown) show a clear difference between the samples from the inner part and from the surface layer of the accumulated snow. All 16 snow samples from the May period taken from the inner part of the snow accumulation have an isotopic composition closer to that of the initial snow compared to the samples from the surface layer of the accumulated snow. In $\mathrm{6}/\mathrm{16}$ experiments, the inner accumulation samples had the same (±uncertainty) isotopic composition for all three isotope variables as the initial snow.

3.2.2 The co-evolution of atmospheric and snow variables

Atmospheric variables

The changes we observe in the PPPs and the isotopic composition of the airborne snow are a result of the interplay between the snow particles and the atmospheric conditions in the wind tunnel. For each experiment, the atmospheric conditions were closely monitored. In the following, we summarise the evolution pattern of all experiments and show the data of experiment no. 9 as a representative example in Fig. 4.

At the start of the experiments, when the conditions inside the wind tunnel were stable, the wind generation was initiated (Fig. 4d), and shortly after, the snow was added to the wind tunnel (Fig. 4, purple shading). Immediately after the snow introduction started, the RH_ice increased sharply and reached saturation with respect to ice within minutes (Fig. 4a). With the start of the wind generator system, the air temperature inside the wind tunnel started to increase as heat was generated by the propeller motor (Fig. 4c, red line). This wind tunnel heating was modulated by the laboratory's cooling system in a way that the cooling phases of the laboratory climate system also slowed down the temperature increase inside the wind tunnel (Fig. 4c, yellow line). The total temperature change for each experiment is given in Table 1. Concurrently with the temperature, the specific humidity inside the wind tunnel evolves as the air reaches (or gets close to) saturation conditions; hence, the saturation vapour pressure defines the humidity level (Fig. 4b). This general evolution of the meteorological variables with the temperature evolution dominating the (saturation) humidity level inside the wind tunnel was observed for all experiments except the warming experiment, in which the air only briefly reaches saturation after ∼ 1 h due to continuous strong heating of the wind tunnel air stream.

The snow-sampling events every 15–30 min can be identified as short-lived dips of 1–2 min in the wind speed (Fig. 4d). In this experiment, the first airborne snow sample (1st $\mathrm{\Delta }{\mathit{\delta }}_{\mathrm{snow}}^{\ast }$) shows a 0.18 ‰ enrichment in δ¹⁸O (Fig. 4f) and a −1.4 ‰ decrease in d-excess (Fig. 4j) whereas the δD value (Fig. 4h) does not change significantly (+0.1 ‰). The maximum observed change (max $\mathrm{\Delta }{\mathit{\delta }}_{\mathrm{snow}}^{\ast })$ for this experiment is a 0.18 ‰ enrichment in δ¹⁸O, a −0.6 ‰ depletion in δD and a −2.0 ‰ decrease in d-excess. The vapour isotopic composition in this experiment showed enrichment in δ¹⁸O (0.81 ‰), δD (10.4 ‰) and an increase in d-excess (3.9 ‰) during the snow introduction. After the snow introduction had ended and until the end of the experiment the δ¹⁸O signal showed a reverse evolution (−0.27 ‰) while both δD (2.9 ‰) and d-excess (5.0 ‰) continued to increase.

Figure 4 The time evolution of all monitored variables in the exemplary experiment no. 9. The meteorological variables inside the wind tunnel are shown in the following panels: (a) relative humidity with respect to ice, (b) CRDS specific humidity readings, (c) wind tunnel and laboratory temperature, and (d) wind speed. The period of snow introduction is shaded in purple. After the snow introduction, the wind tunnel was closed and only opened for brief periods of 1–2 min during airborne snow sampling, which are visible as dips in the wind measurements. The evolution of the 3 min averaged water vapour isotopic composition is shown with shading for the standard error of the mean uncertainty for δ¹⁸O in (e), δD in (g), and d-excess in (i). The isotopic composition of the snow samples is given in (f) for δ¹⁸O, in (h) for δD, and in (j) for d-excess, with the measurement uncertainty as error bars. The snow samples taken from the accumulated snow are shown as square (surface layer/acc_out) and diamond (inside of accumulation/acc_in) markers. The horizontal lines indicate the initial isotopic composition (mean of three samples).

Vapour isotopic composition

The co-evolution of the isotopic composition of the vapour is exemplarily shown in Fig. 4e, g, and i for experiment no. 9 and is summarised for all experiments in Fig. 5 and Table 2, where the observations are divided into changes during and post snow introduction.

During the snow introduction (Fig. 5 a, b, c), most of the 16 experiments in May 2023 showed a strong enrichment in the vapour phase of δ¹⁸O ($\mathrm{13}/\mathrm{16}$) and δD ($\mathrm{16}/\mathrm{16}$) and an increase ($\mathrm{6}/\mathrm{16}$) in d-excess. Two experiments showed a decrease in d-excess. The remaining experiments showed no prominent change during snow introduction in vapour δ¹⁸O ($\mathrm{3}/\mathrm{16}$) or d-excess ($\mathrm{8}/\mathrm{16}$).

Post snow introduction (Fig. 5d, e, f), the evolution in the 14 constant temperature experiments can be categorised into the following:

• i.

  continued increase in δ¹⁸O ($\mathrm{4}/\mathrm{14}$), δD ($\mathrm{13}/\mathrm{14}$), and d-excess ($\mathrm{11}/\mathrm{14}$);

• ii.

  reversed evolution, i.e. strong depletion in δ¹⁸O ($\mathrm{3}/\mathrm{14}$) and δD ($\mathrm{1}/\mathrm{14}$);

• iii.

  low variability in δ¹⁸O ($\mathrm{7}/\mathrm{14}$) and d-excess ($\mathrm{3}/\mathrm{14}$).

In the warming experiment, δ¹⁸O, δD, and d-excess all showed a pronounced increase after the snow introduction had ended until the end of the experiment. Finally, the cooling experiment showed a variable isotopic evolution with a sequence of pronounced phases of depletion and enrichment and another depletion phase in both δ¹⁸O and δD, resulting in net depletion, whereas the d-excess stayed relatively constant (|Δd-excess|<2.7 ‰) during the entire experiment.

Figure 5 The co-evolution of changes in vapour and snow isotopes in all experiments with vapour isotope observations. The change in the snow isotopes (x axes) is plotted against the changes in the vapour isotopes (y axes). The changes in vapour isotopes are calculated from the 3 min averaged data. The upper row (a, b, c) shows changes in vapour isotopes during the snow introduction plotted against the observed change in the snow isotopes of the first airborne sample (see Fig. 3). The lower row (d, e, f) shows the subsequent change in vapour isotopes following snow introduction until the end of the experiment plotted against the maximum observed change in the snow isotopic composition. Shaded areas represent low variability in vapour isotopic composition (|Δδ¹⁸O|<0.3 ‰, |ΔδD|<1.5 ‰, |Δd-excess|<2.7 ‰). Note that the colour code represents the average air temperature during the experiments and allows the comparison between the upper and lower rows and the results shown in Fig. 3. The cooling (black cross) and warming (red cross) experiments are exempt.

4 Discussion

4.1 Airborne snow metamorphism

4.1.1 SSA

All experiments with µCT measurements show a rapidly decreasing SSA with transport time that can be described by the typical exponential decay function (Cabanes et al., 2003). This phenomenon is well known inside the snowpack where the SSA decreases over time in isothermal conditions and under temperature gradient conditions due to snow metamorphism (Harris Stuart et al., 2023; Schleef et al., 2014b), but it seems to happen at an increased rate of 1–2 orders of magnitude during airborne transport. The decrease in SSA describes a change in the particle morphology, and the rapidness with which the SSA decreases suggests a rounding and growing of particles through abrasion or metamorphic mechanisms since large spheres have the lowest SSA values (Legagneux et al., 2002). These results are in line with the changes observed in the sphere size distribution analysis that are discussed below. Note here that simple mechanic fragmentation of snow particles alone does conceptually lead to an increase in sample SSA and thus cannot explain the decrease in SSA. The initial SSA values of the experiments vary slightly, which is due to the storage time between snow production and the experiment start (Table 1). Although the snow was stored at −20 °C to slow down snow metamorphism prior to the experiment start, the older the snow was when introduced into the wind tunnel, the lower the initial SSA was. We find no direct dependency of SSA decay rate on temperature, which is in line with wind tunnel experiments studying the metamorphism rate of accumulated snow under wind influence for temperatures <−5 °C (Walter et al., 2024) but is different from snowpack metamorphism under no-wind conditions (Harris Stuart et al., 2023; Schleef et al., 2014b). The driving processes for snow metamorphism are vapour pressure gradients, which are largely governed by the (absolute) temperature regime in a stationary snowpack (Kaempfer and Schneebeli, 2007; Taillandier et al., 2007). In airborne snow metamorphism, however, other processes such as turbulent mixing of air and vapour, particle–air temperature gradients (see Sect. 4.2), and variability in saturation conditions might dominate the vapour pressure variability and therefore mask the simple absolute temperature dependency that is expected in stationary snowpack metamorphism.

4.1.2 Particle size and shape

The evolution of the fitted sphere size distribution with transport time shows a general broadening and increase in skewness, with a shift of the dominant sphere size to larger diameters. As the algorithm finds the maximum sphere size it can fit at any given location inside the ice matrix of the sample (Hildebrand and Rüegsegger, 1997), this can be interpreted as larger particles emerging and a rounding of the individual particles with transport time. These observations are somewhat surprising since aeolian transport of snow is associated with fragmentation and abrasion of particles and thus a decrease in particle sizes as snow crystals break apart due to mechanical stresses (Comola et al., 2017; Sato et al., 2008). Notably, the evolution of sphere size distribution resembles the evolution observed during isothermal metamorphism (Flin et al., 2004; Legagneux and Domine, 2005). We note here that care must be taken when interpreting fitted sphere sizes as particle sizes (Löwe et al., 2011), and we mention the difficulty in the assessment of particle “size” in snow in general (Fierz et al., 2009). In fact, our observations of growing sphere sizes do not contradict the initial snow particle fragmentation mechanism through impact as pure fragmentation does not significantly influence fitted sphere size distributions within a snow sample. On the other hand, abrasion producing very fine snow dust could increase the occurrence of very small fitted sphere sizes. However, we see a disappearance of the smallest fitted sphere sizes and an appearance of and increase in bigger sphere sizes of >100 µm. One mechanical mechanism that could explain the shift to larger sphere sizes could be the collision and subsequent aggregation of snow particles, as known from ice cloud observations (Pruppacher and Klett, 1997). In ice clouds, this process of coalescing snow particles is an important ingredient in the growth of ice crystals and the formation of snow aggregates (Lo and Passarelli, 1982; Mitchell, 1988).

Another plausible mechanism for the growing of particles is airborne snow metamorphism with continuous sublimation from and vapour deposition on the suspended particles (Walter et al., 2024). The (close-to-) saturation conditions inside the wind tunnel a few minutes subsequent to the snow introduction are favourable conditions for this chain of thermodynamic processes. The growing and rounding of the particles may thus be explained by airborne metamorphic growth, with vapour being preferentially sublimated from convex sub-grain boundaries and entire small grains and preferentially deposited in concavities or on larger grains with lower curvature (Kelvin effect) in analogy to isothermal metamorphism in a snowpack (Colbeck, 1998, 2001; Wakai et al., 2005). Vapour deposition on suspended snow particles has been previously modelled in large-eddy simulations (LESs) of wind-blown snow (Sharma et al., 2018; Sigmund et al., 2024) but could not be validated against observations so far. As was shown by Sharma et al. (2018), these deposition-enabling conditions can exist even under slightly under-saturated surrounding conditions (usually generally favouring sublimation) if the particle temperature is lower than the surrounding air temperature. Traditionally, a temperature equilibrium between the particle and surrounding air is assumed following the well-established and frequently used sublimation rate estimation by Thorpe and Mason (1966). This assumption does not allow for the deposition of vapour on the suspended particle. However, as the temperature in the wind tunnel was not constant, it is possible that the particle temperature was lagging the air temperature, preventing the establishment of a temperature equilibrium between the air and the snow particle. The accumulated snow on the ground and the walls of the wind tunnel could further act as a cooling source in the wind tunnel system that could especially cool saltation particles. Sharma et al. (2018) found a direct dependency between particle diameter and temperature equilibration lag time. This would suggest that, while small particle fragments (from e.g. abrasion) might equilibrate faster and consequently sublimate entirely, bigger and colder particles might attract vapour deposition and thus trigger a positive feedback. Indeed, we find that the smallest sphere sizes detected initially disappear with time, which can be interpreted as the full sublimation of the smallest particle fragments. These explanations combined can resolve the apparently conflicting observations of growing particles due to vapour deposition while we see a net humidity increase from snow sublimation throughout the experiments.

4.2 The isotopic signature of airborne snow metamorphism

We have observed changes in the isotopic composition of vapour and, notably, also in the snow during the experiments. The changes in the snow isotopic composition support the theory of airborne snow metamorphism as a driver for PPP changes since involved phase changes lead to changes in the isotope signals. Mechanical mechanisms such as abrasion, clustering, and aggregation of particles, however, would not lead to changes in the bulk isotopic composition of the snow sample. Investigating the amplitude and sign of the changes in the snow isotopic signal allows further deduction about the processes at play.

4.2.1 Sublimation-controlled metamorphism immediately after snow introduction and during warming experiment

The near-universal signal of enrichment in snow δ¹⁸O and δD while d-excess decreases in the first airborne snow sample suggests sublimation-dominated metamorphic mechanisms affecting the snow and its PPPs during and immediately after snow introduction. This is in line with the steep humidity and RH increase to (close-to-) saturation conditions in this period. Similarly, the continuous positive trend in snow δ¹⁸O and δD and the negative trend in d-excess during the persistent sub-saturated conditions in the warming experiment suggests isotopic fractionation during net sublimation conditions. The existence of fractionation during sublimation has been controversially discussed as the slow process of self-diffusion in ice and usually does not allow mixing on timescales fast enough to incorporate the sublimation signal into the remaining solid matrix (Friedman et al., 1991; Hu et al., 2022). Thus, it was assumed that sublimation happens according to a layer-by-layer removal of the water molecules, which would not alter the isotopic composition of the remaining ice. Following this reasoning, current isotope-enabled climate models do not account for fractionation during snow sublimation (Cauquoin and Werner, 2021; Dütsch et al., 2023; Moore et al., 2016; Pfahl et al., 2012). However, ever-growing evidence documents an enrichment in δ¹⁸O and δD in snow isotopic composition during and after sublimation of surface snow (Hughes et al., 2021; Sokratov and Golubev, 2009; Stichler et al., 2001) analogously to the results of Δδ¹⁸O and ΔδD in the first airborne snow samples and the trend analysis in the warming experiment. The concurrent decrease in snow d-excess is also in agreement with previous laboratory studies that likewise found a decrease in d-excess under sublimation conditions (Hughes et al., 2021). It is an indicator of the existence of non-equilibrium, diffusive processes, which aligns with the strong sublimation signal in the humidity observations. Thus, the wind tunnel results presented in this study further add to previous studies with δ¹⁸O and δD enrichment and d-excess decrease during sublimation-favourable conditions.

These sublimation-dominated conditions in the wind tunnel at the beginning of the experiments can explain the PPP changes observed in the first airborne samples. The shift in the smallest sphere sizes is largest between the initial snow and the first airborne snow samples. Highly dendritic parts of the airborne particles either first get abraded and fragmented or sublimate directly from the airborne dendritic particles, which explains the disappearance of the smallest fitted sphere sizes and the decrease in SSA. Once (close-to-) saturation conditions are reached, the conditions shift towards deposition-controlled metamorphism.

4.2.2 Deposition-controlled metamorphism during aeolian transport

We hypothesise that the further decrease in SSA and the gradual shift towards bigger sphere sizes is due to vapour deposition on the airborne particles that consequently grow with aeolian transport time. This hypothesis is supported by the observed maximum changes in snow isotopic composition that can be explained by assuming vapour deposition on suspended particles. Following this assumption, the vapour that gets deposited on the suspended particles controls the sign and magnitude of the change observed in the airborne snow samples. Indeed, the changes can be explained by analysing the isotopic composition of the solid that is formed (deposited) from vapour inside the wind tunnel and the initial snow that is introduced. Figure 6b shows this analysis for δD assuming equilibrium fractionation. The solid disequilibrium analysis (Sect. 1) shows that, generally, experiments in which the combination of water vapour and temperature inside the wind tunnel produces an equilibrium deposit that is enriched compared to the initial snow (positive disequilibrium) show a positive change in δD. If, however, the formed equilibrium deposit is depleted compared to the initial snow (negative disequilibrium), the change in δD is negative. A linear fit to this data is significant (p=0.004), with R²=0.46. It is worth mentioning that a temperature uncertainty of +1 K can produce a fit with a zero origin. In Fig. 6a, the same analysis is shown for δ¹⁸O. As δ¹⁸O is more sensitive to non-equilibrium conditions than δD (Majoube, 1970; Merlivat and Nief, 1967), it is not surprising that the disequilibrium analysis of formed equilibrium deposits cannot explain the observations. Instead, the changes in δ¹⁸O can be explained by allowing for a colder particle temperature compared to the surrounding air. As discussed earlier, this condition of a particle–air temperature gradient that leads to supersaturation with respect to the particle surface is a requirement for the growth of particles through vapour deposition on such short timescales, and the changes in δ¹⁸O support this theory. Isotopic fractionation under supersaturation conditions can produce solid deposits with very depleted δ¹⁸O values compared to pure equilibrium fractionation (Jouzel and Merlivat, 1984). The resulting δ¹⁸O value of the deposit is dependent on the supersaturation condition (set by the particle–air temperature difference) that strongly influences the non-equilibrium fractionation of δ¹⁸O. The changes in snow δ¹⁸O during the wind tunnel experiments are thus dependent on the combination of initial snow and vapour isotopic composition and the particle–air temperature difference that defines the supersaturation. Everything except the particle–air temperature difference is measured in the wind tunnel experiments, which allows us to estimate the temperature difference needed to observe negative changes in snow δ¹⁸O values. Following Jouzel and Merlivat (1984), we calculate the minimum supersaturation conditions necessary for negative changes in snow δ¹⁸O. From this, we deduce the associated temperature difference between air and particles using the saturation vapour pressure functions of Goff and Gratch (1946) (see colour code in Fig. 6a). As can be seen, the four experiments in which a negative change in snow δ¹⁸O was observed during the experiments are indeed cases in which reasonable temperature differences of 0.5–1.0 K are enough to result in the negative δ¹⁸O changes that were observed. Generally, the calculated required temperature difference shows a strong correlation with the maximum observed changes in snow δ¹⁸O, with R²=0.47 and p=0.004 (not shown).

In summary, the observed changes in snow isotopic composition support the theory of growing particles due to vapour deposition and the existence of the prerequisite of particles that are colder than the surrounding air.

Figure 6 Changes in snow isotope signal explained. The y axes show the maximum isotopic change observed between the airborne snow sample and the initial-snow sample for δ¹⁸O in (a) and for δD in (b) (same as Fig. 3d and e). The x axes show the average solid disequilibrium in the period after the snow had been added until the maximum change in snow isotopic composition was observed. The colour code in (b) shows the average temperature of the experiments. In (a), the colour code shows the minimum particle–air temperature difference required to produce supersaturation conditions that lead to a negative δ¹⁸O change in the snow. Linear least-square fits to the data and the corresponding fitting statistics are given.

4.2.3 Fractionation during sublimation?

As outlined above, this study finds isotopic enrichment in the first airborne snow samples under net sublimation conditions, supporting previous snow-focused studies that present evidence for fractionating sublimation (Hughes et al., 2021; Sokratov and Golubev, 2009; Stichler et al., 2001). A recent study by Wahl et al. (2021) further substantiated the existence of isotopic fractionation during net sublimation from a snow surface by measuring the isotopic composition of the near-surface sublimation flux, which was shown to be significantly depleted compared to the surface snow isotopic composition. Besides the change in the snow isotope signal, the wind tunnel data provide additional supporting evidence through the vapour isotope signal. That is, the changes observed in the vapour isotopic composition during snow introduction are better described by assuming a mixing between already present and newly produced water vapour originating from fractionating sublimation rather than from non-fractionating sublimation. Specifically, the initial disequilibrium prior to snow introduction (i.e. vapour produced by fractionating sublimation) is a better predictor for the changes in vapour isotopic composition during the snow introduction (R²(δ¹⁸O) = 0.76, R²(δD) = 0.68) than the mere difference in vapour and snow isotopic composition (i.e. vapour produced by non-fractionating sublimation) (R²(δ¹⁸O) = 0.47, R²(δD) = 0.28). For corresponding plots, see SM Fig. S3. However, we want to point out that our observations of an isotopic-deposition signal during general net sublimation (specific humidity increase) suggests an alternative explanation for the observed fractionation during net sublimation. Snow metamorphism happens on the micro-scale, driven by net vapour–pressure gradients where “metamorphic sublimation” and “metamorphic deposition” can co-exist within the snowpack while being spatially separated on the micro-scale, e.g. within one snow crystal. This process happens continuously in the snowpack (Colbeck, 1982) and, evidently, also during airborne transport, with periods of varying dominance of one process over the other. This concept of airborne snow metamorphism is schematically visualised in Fig. 7. Thus, observations of isotopic fractionation during net sublimation do not require that metamorphic sublimation is indeed fractionating if simultaneous metamorphic deposition is. The isotopic enrichment in the solid matrix could thus stem from enriched deposits rather than a surface becoming enriched by sublimation. As in our experiment, this would require that the deposit has a lower d-excess than the snow particle on which it gets deposited, leading to a d-excess decrease rather than the theoretical d-excess increase during vapour deposition (Casado et al., 2016). Ultimately, both d-excess increases and decreases in a snowpack can be associated with metamorphic deposition depending on the disequilibrium between snow and vapour. The macro-scale humidity gradient, e.g. the gradient between the snow surface and the near-surface atmosphere, defines the direction of the net humidity flux, but isotopic fractionation is defined through the micro-scale metamorphic processes. These considerations are needed to find a physically sound explanation for the evidently fractionating net sublimation. Clearly, for realistic isotope-enabled climate model simulations, an implementation of a bulk fractionation factor for net surface sublimation is needed.

Figure 7 Conceptual scheme of airborne snow metamorphism. The scheme shows that micro-scale processes of metamorphic sublimation and metamorphic deposition happen continuously and simultaneously even under net sublimation conditions. We hypothesise that metamorphic vapour deposition on the largest airborne particles due to temperature differences between the surrounding air and particle surface, together with preferential sublimation of convex sub-grain boundaries and preferential deposition in concavities, causes a decrease in SSA and particle growth.

4.3 Field application and experimental limitations

Airborne snow metamorphism as observed in the ring wind tunnel experiments is expected to be dominant in the saltation layer of wind-blown snow. The (close-to-) saturation conditions inside the ring wind tunnel better resemble saltation layer conditions than suspension layer conditions on average, but Sigmund et al. (2021) have shown that downward moisture fluxes can also occur in suspension clouds. Saturated saltation layers (<1 m height) have been frequently observed in coastal Antarctic areas, whereas saturation in blowing-snow clouds of heights >5 m requires very high wind speeds and high particle concentrations (Amory and Kittel, 2019). A deepening of the saturated layer is favoured in stable atmospheric boundary layers with strong near-surface temperature gradients (Mann et al., 2000). Thus, we expect airborne snow metamorphism to reach larger heights during stable boundary layers, which are frequently observed during polar-night conditions. Additionally, the particle–air temperature gradient needed to facilitate vapour deposition and particle growth is predominantly found in the saltation layer, whereas the thermal-equilibrium assumption seems to adequately describe the conditions in the suspension layer (Sharma et al., 2018). Indeed, the particle–air temperature difference simulated in LES close to the cold snow surface is of comparable magnitude (∼ 0.5 K; Sigmund et al., 2024) to the results obtained from the δ¹⁸O isotope analyses. The sustained deposition conditions in our experiments might also be favoured by the continued warming of the wind tunnel air, with the particle temperatures lagging. The motor heating is an experimental limitation that hinders the direct analogy to real-world conditions where the energy balance inside a blowing-snow cloud is largely defined by the sublimation process itself, with sublimation leading to a cooling which essentially self-inhibits ongoing sublimation (Amory and Kittel, 2019; Sigmund et al., 2021).

An additional experimental limitation is the interaction of airborne particles with the physical boundaries of the wind tunnel, i.e. the walls and floor, which results in additional collision impacts compared to in nature or in straight wind tunnels. Furthermore, the flow field inside the wind tunnel does not fully develop into a logarithmic boundary layer flow, and centrifugal forces act on the particles in the curved section (Sommer et al., 2017). However, as Walter et al. (2024) outline for the same wind tunnel, the individual particle transport characteristics of the impact angle and velocities, both on the wind tunnel walls and at the snow surface, are comparable to particle transport characteristics in the saltation layer of a fully developed boundary flow. In summary, the sheer setup of the wind tunnel construct prevents an unequivocal comparison to real-world conditions. Yet, the experiments have provided a means to discover a physical process that is not discernible from fixed-location in situ measurements of wind-blown snow without knowledge of the provenance and prior evolution of the observed airborne particles. We point out that the absolute values of available field observations of particle shape and size distributions do not conflict with the findings of this study but might have implicitly entailed airborne snow metamorphism markings. The largest particles that are observed in the experiments contain sphere sizes of 200–300 µm. These sphere sizes are well within in situ observed particle size distributions (Nishimura et al., 2014; Nishimura and Nemoto, 2005) and are significantly smaller than fresh snow particles (Comola et al., 2017; Sato et al., 2008), therefore not contradicting snowflake fragmentation. Interestingly, the most abundant sphere sizes in the last airborne snow samples have diameters of 60–70 µm, which would be considered to be small particle sizes in wind-blown particle size distributions. However, care must be taken when comparing particle size estimates from snow particle counter (Nemoto and Nishumura, 2001) and sphere-fitting routines (Hildebrand and Rüegsegger, 1997). Which metric is superior depends on the application and the question that is being asked. Highly dendritic fresh snowflakes can be large in terms of the particle size but can have very high SSA values, which is better represented through the small sphere size metric.

4.4 Implications of airborne snow metamorphism

4.4.1 Influence on snow surface and snowpack characteristics

The wind tunnel experiments have shown that airborne particles can grow thermodynamically under certain conditions. This process is not accounted for in wind-blown snow models (Melo et al., 2022), nor is it accounted for in calculations of mass removal by snow sublimation during blowing-snow events (Palm et al., 2017). Unfortunately, we have no means to estimate the mass of redeposited vapour in our experiments. Estimates of the impact on mass balance and sublimation flux calculations are beyond the scope of this study. Yet, the changes in airborne PPPs can be set in the context of known changes in snowpack characteristics under wind influence.

Wind has been shown to increase the density of snow (Lehning et al., 2002; Liston et al., 2007; Walter et al., 2024). Crucial for this is the presence of wind-transported snow (Sommer et al., 2017). The surface snow density increase is associated with grain settling and packing, which are more efficient for smaller and rounder particles (Parteli et al., 2014; The Firn Symposium team et al., 2024). The combination of mechanical abrasion and fragmentation (Comola et al., 2017), together with metamorphic rounding, could thus increase the densification rate under wind influence. Likewise, SSA has been shown to decrease under wind influence (Legagneux et al., 2002; Walter et al., 2024). We suggest that this observed wind–SSA relationship is largely due to airborne snow metamorphism. These modifications of SSA affect the radiative characteristics of the resulting snowpack, with SSA decay leading to a reduction in snow albedo (Flanner and Zender, 2006). Consequentially, airborne snow metamorphism could reduce the albedo of the ultimately accumulated snow.

4.4.2 The isotopic signature of wind-blown snow

Snowpack metamorphism has been shown to influence the isotopic composition of snow in polar areas depending on the temperature gradient in the snowpack (Casado et al., 2021; Harris Stuart et al., 2023). Laboratory experiments have further shed light on the interplay between snow and vapour isotopic composition during snow metamorphism (Ebner et al., 2017). The experiments revealed that the isotopic signature of the wind-blown snow event in the wind tunnel was dependent on the combination of initial vapour and snow isotopic composition and the resulting disequilibrium. As we produced snow from tap water, this disequilibrium was, at times, asymmetric between the two isotope species and was more pronounced than generally assumed in nature. As the disequilibrium determines the isotopic evolution, it is not possible to unambiguously predict the expected changes in the snow isotope signal under wind influence without considering the water vapour isotope variability imposed by synoptic-scale atmospheric transport (Aemisegger et al., 2022; Bagheri Dastgerdi et al., 2021). However, the results suggest that a strong d-excess decrease can be linked to airborne metamorphism. This should be kept in mind when observations of snow d-excess values are used as hydrological tracers. Further field studies in windy and dry locations, such as the katabatic wind zones on the Antarctic Ice Sheet, could support this idea.

Another interesting implication of the wind tunnel findings is the interpretation of the snow isotopic composition as a climate signal. Stable water isotope records obtained from ice cores in the polar areas are valuable climate proxy records that allow the reconstruction of paleo-climate and ice sheet dynamics in these remote areas (Jouzel et al., 2007; Stenni et al., 2017; Vinther et al., 2009). The foundation for these paleo-climate interpretations of isotope records is the observation of a linear relationship between annual mean temperatures and the corresponding average precipitation isotopic composition (Dansgaard, 1964). However, in low-accumulation, dry regions of the ice sheets, from which the oldest ice cores are obtained, evidence accumulates for the influence of post-depositional processes over the snow isotopic composition and thus over the climate signal that is stored in the snow (Dietrich et al., 2024; Ollivier et al., 2024; Town et al., 2024; Wahl et al., 2022; Zuhr et al., 2023). This study shows that wind transport is one more post-depositional process that has the potential to add a local climate signal to the snow isotopic composition. The magnitude of observed snow isotopic change in the experiments due to airborne snow metamorphism is within the range of observed daily snow surface isotope variability on clear-sky and calm days (Wahl et al., 2022; Zuhr et al., 2024), yet this evolves within 1–2 h. Thus, modelling studies will be needed to estimate the impact of the combined post-depositional processes on the original precipitation isotopic composition that is also proclaimed to carry information about the moisture source conditions (Aemisegger et al., 2022; Uemura et al., 2012). Our results suggest that locations with frequent strong winds might show strong potential for an overprinting of the original d-excess source region signal and could support the interpretation of studies with measured low d-excess values (Hu et al., 2022).

5 Conclusions

We used ring wind tunnel experiments to study the governing micro-scale processes of aeolian snow transport that are difficult to observe in the field. Micro-CT measurements of airborne snow samples documented a growing and rounding of suspended particles and a consequential decrease in SSA with transport time. Stable water isotope analysis of airborne snow and water vapour suggests airborne snow metamorphism rather than mechanical mechanisms as a driver for the observed particle evolution. Snow sublimation dominated the metamorphic process during sub-saturation conditions, while vapour deposition dominated in (close-to-) saturated conditions comparable to saltation layer conditions during wind-blown snow events. Supporting previous modelling studies, the wind tunnel results further indicated a particle–air temperature difference, particularly for larger particles that favour the deposition of vapour and thus particle growth, while smaller particles get sublimated entirely. This process of airborne snow metamorphism is difficult to observe in situ with point measurements of wind-blown snow, where the origin and preceding history of the airborne particles are unknown. Thus, it is a hitherto unobserved process, yet our conclusions do not contradict available in situ measurements of airborne particle distributions during blowing-snow events. Instead, our findings refine the process understanding of wind-blown snow, with possible implications for mass and energy balance and moisture budget calculations during blowing- and drifting-snow events.

In this study, the isotopic fingerprint of airborne snow metamorphism was used to pinpoint its existence. The change in the snow isotope signal that we attribute to airborne snow metamorphism was dependent on the vapour–snow disequilibrium and ranged from −0.88 ‰ to +1.47 ‰ in δ¹⁸O, −5.7 ‰ to +5.7 ‰ in δD, and −6.1 ‰ to +1.3 ‰ in d-excess. As this isotopic signature of airborne snow metamorphism partly overprints the original snow isotopic signal, it offers the potential to act as tracer for aeolian snow transport. Furthermore, it questions the conservation of moisture source region signals in the snow isotopic composition in windy areas with frequent wind-blown snow events. Future field-based studies will be valuable to support our laboratory-based observations of a link between airborne snow transport, particle physical properties, and stable water isotopes.