On the importance of sublimation to an alpine snow mass balance in the Canadian Rocky Mountains

A modelling study was undertaken to evaluate the contribution of sublimation to an alpine snow mass balance in the Canadian Rocky Mountains. Snow redistribution and sublimation by wind, snowpack sublimation and snowmelt were simulated for two winters over an alpine ridge transect located in the Canada Rocky Mountains. The resulting snowcover regimes were compared to those from manual snow surveys. Simulations were performed using physically based blowing snow (PBSM) and snowpack ablation (SNOBAL) models. A hydrological response unit (HRU)based spatial discretization was used rather than a more computationally expensive fully-distributed one. The HRUs were set up to follow an aerodynamic sequence, whereby eroded snow was transported from windswept, upwind HRUs to drift accumulating, downwind HRUs. That snow redistribution by wind can be adequately simulated in computationally efficient HRUs over this ridge has important implications for representing snow transport in large-scale hydrology models and land surface schemes. Alpine snow sublimation losses, in particular blowing snow sublimation losses, were significant. Snow mass losses to sublimation as a percentage of cumulative snowfall were estimated to be 20–32% with the blowing snow sublimation loss amounting to 17–19% of cumulative snowfall. This estimate is considered to be a conservative estimate of the blowing snow sublimation loss in the Canadian Rocky Mountains because the study transect is located in the low alpine zone where the topography is Correspondence to: J. W. Pomeroy (john.pomeroy@usask.ca) more moderate than the high alpine zone and windflow separation was not observed. An examination of the suitability of PBSM’s sublimation estimates in this environment and of the importance of estimating blowing snow sublimation on the simulated snow accumulation regime was conducted by omitting sublimation calculations. Snow accumulation in HRUs was overestimated by 30% when neglecting blowing snow sublimation calculations.


Introduction
Snowpack depth and density in the alpine zone of high mountains exert a strong control on the magnitude, timing and duration of snowmelt as well as directly influencing alpine ecology and avalanche formation (Jones et al., 2001).Snowcover increases the surface albedo and provides a colder surface to interact with the atmosphere compared to snow-free zones.Thus, there are marked differences in energy and moisture fluxes over snow-covered and snow-free surfaces, which have implications for evapotranspiration, permafrost, and glaciers.Snowpack and snowcover characteristics in alpine zones are strongly influenced by wind through the action of wind in entraining, transporting and sublimating snow (Dyunin and Kotlyakov, 1980;Föhn, 1980;Schmidt et al., 1984;Meister, 1989;Pomeroy, 1991).
Snow sublimation has been shown to be an important part of alpine water balances (e.g.Kattelmann and Elder, 1991;Marks et al., 1992;Strasser et al., 2008).Since snowmelt has an extremely important role in regulating annual runoff, these winter sublimation losses can have an important effect 1402 M. K. MacDonald et al.: The importance of sublimation to an alpine snow mass balance on the hydrology of cold regions.Sublimation of snow is the conversion of snow particles from solid to vapour, bypassing the liquid water phase.Snow sublimation requires an energy input of 2.84 MJ kg −1 (at 0 • C and standard atmospheric pressure) and the vapour pressure at the snow surface must be greater than that of the ambient air.Snow sublimation rate calculations are based on Thorpe and Mason's (1966) relationships for the sublimation rate of an ice sphere based on a balance of atmosphere-particle heat transfer and latent heat due to sublimation at the particle surface.Mass loss due to sublimation occurs, in order of decreasing efficiency, from wind transported snow, intercepted snow and from the snowpack.Strasser et al. (2008) showed that cumulative sublimation from snowpack, intercepted snow in forest canopies and blowing snow amounted to 24% of snowfall over a one-year period in a mountainous region in southeast Germany.Over 70% of annual snowfall sublimated during wind transport at certain summit crest and ridge locations.
Wind transport of snow is a common phenomenon across high altitude and latitude cold regions that can significantly affect snowcover distribution patterns both during accumulation-and ablation-dominated periods.Snow transport involves the horizontal redistribution of snow.Surface snow is eroded and transported via saltation (Schmidt, 1986;Pomeroy and Gray, 1990) and suspension (Budd et al., 1966;Pomeroy and Male, 1992) from flat surfaces, hilltops, windward slopes and sparsely vegetated surfaces to topographic depressions, leeward slopes and more densely vegetated surfaces (Pomeroy et al., 1993;Liston and Sturm, 1998).The occurrence of blowing snow is dependent upon the local boundary layer meteorology and physical characteristics of the snow surface.Blowing snow can proceed when the surface wind speed exceeds a threshold wind speed dictated by the snow cohesive bond forces, which depend on the thermal and settling histories of the snowpack (Li and Pomeroy, 1997a).Snow particles transported by wind are well ventilated and undergo sublimation in the presence of an unsaturated atmosphere (Dyunin, 1959;Schmidt, 1972Schmidt, , 1986).The blowing snow particle sublimation rate is controlled by atmospheric turbulence, temperature, humidity, incoming radiation and particle size.Blowing snow sublimation losses of 15 to 41% of annual snowfall have been estimated for the Canadian Prairies (Pomeroy and Gray, 1995), 28% of annual snowfall over western Canadian Arctic tundra (Pomeroy et al., 1997), 18-25% of winter precipitation over the Alaskan arctic (Liston and Sturm, 2002), up to 20% of the annual precipitation over certain areas of the Antarctic ice sheet (Bintanja, 1998) and only 3.7% of snowfall during an austral winter period at Halley Station in Antarctica (Mann et al., 2000).There are two general types of model representations of blowing snow sublimation that provide considerably different estimates of blowing snow sublimations quantities.One type calculates the sublimation of blowing snow particles using observed temperature and humidity forcing, an expression for ventilation velocity as a function of friction velocity and an expression for the vertical gradient in undersaturation during blowing snow events (Pomeroy et al., 1993;Liston and Sturm, 1998).Other models allow a vertical humidity profile to develop by calculating a thermodynamic feedback (i.e. a negative feedback in the form of decreased water vapour deficit and temperature) during the sublimation of blowing snow (Bintanja, 1998;Bintanja and Reijmer, 2001;Déry et al., 1998;Déry and Yau, 1999;Déry and Yau, 2001).This negative feedback in turn reduces calculated sublimation rates.There is still debate over this model component as some argue that atmospheric mixing during blowing snow events entrains sufficient dry air so as to diminish the effect of the negative feedback (Pomeroy and Li, 2000;Bintanja, 2001).The snow transport models presented by Doorschot et al. (2001) and Lehning et al. (2008) neglect the sublimation of blowing snow particles.
The interception of snowfall by forest canopies is hydrologically important as intercepted snow is susceptible to greater sublimation rates than snow on the ground.Though snow interception by forest canopies and subsequent sublimation is not a component of alpine snow mass balances, it can be an important component of treeline snow mass balances adjacent to alpine zones.The sublimation of intercepted snow can be substantial as the ratio of snow surface area to mass is large, allowing for the turbulent transfer of heat to intercepted snow to be sufficient for sublimation to occur.Intercepted snow sublimation losses have been estimated to be 20-30% of snowfall at a subalpine forest site in the continental United States (Montesi et al., 2004), only 4.5-5.2% of snowfall in northern Idaho (Satterlund and Haupt, 1970), and at a rate of 0.71 mm day −1 at a subalpine forest site in the Colorado Rocky Mountains (Molotch et al., 2007).
The sublimation (or condensation) of snow on the ground is ubiquitous where seasonal snowcover is present, though the sublimation rate is small relative to that for blowing snow and for intercepted snow due to a smaller ratio of surface area to mass and a lower exposure to atmospheric ventilation.Surface snow sublimation losses during winter periods have been estimated to be 0.36 mm d −1 in continental Sweden (Bengtsson, 1980), 1-2 mm d −1 in the eastern Canadian Rocky Mountains during Chinook events (Golding, 1978), and 15% of snowfall at an alpine site in the Colorado Front Range (Hood et al., 1999).
Hydrological and atmospheric models require some description of blowing snow redistribution and sublimation that is suitable for complex terrain for application to cold regions (Dornes et al., 2008).The large scale application of these models in mountain and polar environments precludes a finely distributed approach such as those employed for small basins (e.g.Liston and Sturm, 1998;Essery et al., 1999) and some form of landscape aggregation is necessary.Given suitable single column estimates of snow transport and sublimation, a particular difficulty lies in discretizing blowing snow erosion and deposition zones.The accurate estimation of enhanced snow accumulation in drift zones is highly sensitive to both the snow mass available for deposition (controlled by the area of the blowing snow erosion zone) and the area over which the snow mass is deposited.Approaches haven been developed in the context of avalanche forecasting (Durand et al., 2001;Lehning and Fierz, 2008).Durand et al. (2001) developed an operational parameterization for snow redistribution over mountainous terrain, though an evaluation of a spatially distributed application over actual terrain has not been presented.Lehning and Fierz (2008)  The purpose of this study is to evaluate the contribution of sublimation to the winter alpine snow mass balance in the Canadian Rocky Mountains, with an emphasis on its impact on the end-of-winter snow accumulation that largely controls snowmelt runoff and the spring freshet in this headwater region.A particular focus will be on blowing snow since 1) it is the most efficient snow sublimation process, 2) model representations of blowing snow sublimation vary considerably and 3) estimating snow accumulation regimes due to snow redistribution is highly sensitive to the spatial discretization of the model domain.The specific objectives are to identify HRUs that are suitable for simulating snow accumulation and redistribution over alpine topography in mountains with strong winds, and to estimate the winter snow sublimation loss at an alpine site in the Canadian Rockies.The test area is in the Front Ranges of the Canadian Rocky Mountains.which is characterized by sharp topographic gradients and steep slopes, Chinook winds and other strong wind features.

Site description
Fisera Ridge (hereafter FR; 50 • 57 N; 115 • 12 W; 2305-2322 m a.s.l.) is an alpine ridge located within the Marmot Creek Research basin (MCRB; Fig. 1).MCRB is a 9.5 km 2 watershed located in the Rocky Mountain Front Ranges in Alberta, Canada.The general aspect is easterly.The basin is primarily montane with subalpine forest and alpine tundra ridgetops.The basin landcover consists of dense lodgepole pine, mature spruce and subalpine fir forest at lower elevations, larch, shrubs and grasses at and just below the treeline, and talus and bare rocks in the high alpine.MCRB is underlain by glacial and post-glacial deposits ranging from 10 to 30 m depth above bedrock, except at high elevations and along portions of the creek channels (Stevenson, 1967).Seasonally frozen soils are present at higher elevations.Annual precipitation is typically around 900 mm with 60-75% being snow.Climate normals as recorded at the Kananaskis Pocaterra station (ID 3053604; 1610 m a.s.l.) range from a low of −11.1 • C in January to a high of 11.4 • C in August.Temperatures are typically colder at MCRB since it is at a higher elevation.Marmot Creek itself is a tributary of the Kananaskis River and is part of the Bow River system.FR is located above the treeline, where subalpine fir and larch give away to sparse shrubs, exposed soils and grass.The highest elevation of FR has a western boundary with an elevation of approximately 2617 m a.s.l., and decreases in elevation  1).The ridge-top station is located at the top of FR and measures air temperature, relative humidity, incoming shortwave radiation, incoming longwave radiation, wind speed and direction.The north-facing and southeast-facing stations are located on the northern and southern faces of FR, respectively, and measure wind speed and direction.Wind measurements were made at approximately 2 m above the ground surface.Average measured wind speed, temperature and relative humidity over the simulation periods are presented on Table 1.A Geonor T200B accumulating precipitation gauge was installed in a sheltered area near the ridge-top FR station during the fall of 2008.Thus, for 2008/2009, precipitation data from the FR Geonor gauge was used.Elevation-induced differences in precipitation exist between the FR and UC sites.2008/2009 FR precipitation data was correlated with precipitation data from the UC Geonor T200B accumulating precipitation gauge to develop a multiplier (1.18) to extrapolate 2007/2008 UC precipitation data to the FR site.The Geonor precipitation gauge data were corrected for undercatch according to the equation presented by MacDonald and Pomeroy (2007).Atmospheric pressure is measured at the HM station.Standard atmospheric pressure relationships and known site elevations were used to estimate FR atmospheric pressure from the HM observations.Manual snow surveys were performed over FR during 2007/2008 and 2008/2009.The snow survey transect extended 200 m from the NF station over FR, beyond the SF station and into the forested area (Fig. 2).The modelled do-main corresponds to the snow survey transect.Snow depth was measured every 1-3 m and snow density was measured every fifth depth measurement using an ESC30 snow tube when possible.The snowpack was often too shallow to measure on the windswept north-facing slope and too deep (>120 cm) to measure on the south-facing slope with the ESC30 tube.Snow pits were dug when possible at the locations shown on Fig. 2. Snow density was measured in the snow pit to depth by weighing samples obtained using a fixed triangular cutting device (Perla "Swedish Sampler").To calculate mean snow water equivalent (SWE) for an HRU, the mean measured snow density for a particular HRU was applied to each depth measurement in that HRU.
A vegetation survey was conducted along the FR snow survey on 3 July 2008.A shrub count was performed within two 9 m×9 m grids (one on the north-facing slope and one on the south-facing slope).Eight shrubs were counted within the north-facing slope grid and 47 shrubs were counted within the southeast-facing slope grid, yielding 0.1 shrubs m −2 and 0.6 shrubs m −2 on the on the northfacing and southeast-facing slopes, respectively.Twentythree shrub height and width measurements were taken along the snow survey.Shrub measurements were performed manually with a ruler.The shrub height measurement was taken to be the distance from the ground surface to the maximum vertical extent of the shrub.The shrub width measurement was taken to be the maximum girth of the shrub.Average shrub height and width along the transect were 63 cm and 108 cm, respectively.Average shrub height was 51 cm and 82 cm on the north-facing and southeast-facing slopes, respectively.Average shrub width was 107 cm and 111 cm on the north-facing and southeast-facing slopes, respectively.The shrub width measurements included the aggregation of several clumps of shrubs.Photographs taken with a camera equipped with a hemispherical lens were analyzed with GLA software (Frazer et al., 1999) to determine the leaf area index (LAI) for the spruce forest downslope from the southeast-facing station.An average LAI of 0.91 and an average canopy height of 2.3 m was determined.
An airborne light detection and ranging (LiDAR) data collection campaign was deployed over MCRB research during August 2007.High-resolution digital elevation data was obtained.A 10 m DEM of MCRB was created using this highresolution LiDAR data using Golden Software Surfer 8.00.

Models used
A suite of physically based algorithms were used to simulate snow accumulation over FR.The algorithms were combined within the Cold Regional Hydrological Modelling platform (CRHM; Pomeroy et al., 2007).CRHM is an object-oriented hydrological modelling platform developed for Canadian environments (e.g.boreal forest, mountain forests, alpine tundra, muskeg, arctic tundra and prairies).The spatial discretization is in the form of HRUs with a conceptual landscape sequence or water flow cascade.CRHM has a modular structure in that a model is created by selecting from a library of process modules.The following CRHM modules were used in this study: Global and Slope Qsi (radiation calculations with adjustments for aspect, elevation and slope), PBSM (snow transport and sublimation), SNOBAL (snowmelt), Canopy (adjusts shortwave and longwave radiation exchanges beneath needleleaf forest canopies and accounts for canopy effects on water mass balance at the ground surface).
The snow mass balance over a uniform element of a landscape (Fig. 3) is the result of snowfall accumulation, the distribution and divergence of blowing snow fluxes both within and surrounding the element, and sublimation and melt from the snowpack.The terms presented on Fig. 3 are described in the CRHM module subsections.

Global and Slope Qsi
The CRHM Global module calculates theoretical clear-sky direct and diffuse solar radiation.Global calculates the theoretical direct beam component of solar radiation, Q dir , according to Garnier and Ohmura (1970) and the diffuse clearsky radiation component, Q dif , according to List (1968) as (1) where I is the intensity of extraterrestrial radiation, p is the mean zenith path transmissivity of the atmosphere, m is the optical air mass (calculated from Kasten and Young, 1989), δ is the declination of the sun, θ is the latitude, H is the hour angle measured from solar noon positively towards west, A is the slope azimuth measured from the north through east, Z is the slope angle, a w is the radiation absorbed by water vapour (7%), a c is the radiation absorbed by ozone (2%) and Q ext is the extraterrestrial radiation on a horizontal surface at the outer limit of the earth's atmosphere.The Slope Qsi module calculates incident solar radiation on slopes based on the ratio of measured incident shortwave radiation on a level plane and the calculated clear-sky direct and diffuse shortwave radiation on a level plane (from Global).

Prairie blowing snow model
PBSM calculates two-dimensional blowing snow transport and sublimation rates for steady-state conditions over a landscape element using mass and energy balances.PBSM was initially developed for application over the Canadian Prairies, characterized by relatively flat terrain and homogeneous crop cover (e.g.Pomeroy, 1989;Pomeroy et al., 1993).Certain assumptions and parameterizations in PBSM were derived from field observations in the Canadian Prairies and therefore should be applied outside this environment with caution.However, versions have been applied to variable vegetation height (Pomeroy et al., 1991), over alpine tundra (Pomeroy, 1991), arctic tundra (Pomeroy and Li, 2000) and mountainous subarctic terrain (MacDonald et al., 2009) and show an ability to simulate winter snowpack evolution.The model is used here because it has well tested comprehensive snow transport and sublimation calculations and considers the effect of vegetation.Full details of the model are presented by Pomeroy and Gray (1990), Pomeroy and Male (1992), Pomeroy et al. (1993) and Pomeroy and Li (2000).
The snow mass balance over a uniform element of a landscape (e.g. a HRU) is a result of snowfall accumulation and the distribution and divergence of blowing snow fluxes both within and surrounding the element given by where dS/dt is the surface snow accumulation (kg m −2 s −1 ), P is snowfall (kg m −2 s −1 ), p is the probability of blowing snow occurrence within the landscape element, F is the blowing snow transport out of the element (kg m −1 s −1 ) which is the sum of snow transport in the saltation and suspension layers, F salt and F susp , ∇ H is the horizontal divergence, E B (x)dx is the vertically integrated blowing snow sublimation rate (kg m −1 s −1 ) over fetch distance x (m), E is the snowpack sublimation (kg m −2 s −1 ) and M is snowmelt (kg m −2 s −1 ).
Since PBSM is for fully-developed blowing snow conditions, PBSM is restricted to minimum fetch distances of 300 m following measurements by Takeuchi (1980).Blowing snow transport fluxes are the sum of snow transport in the saltation and suspension layers, F salt and F susp (kg m −1 s −1 ), respectively.Saltation of snow must be initiated before snow transport can occur in the suspension layer and blowing snow sublimation can occur.
F salt is calculated by partitioning the atmospheric shear stress into that required to free particles from the snow surface, to that applied to nonerodible roughness elements (exposed vegetation) and to that applied to transport snow particles (Pomeroy and Gray, 1990).Estimation of the threshold condition and shear stress required to free particles from the surface follows the procedure of Li and Pomeroy (1997).The shear stress applied to vegetation is estimated following Pomeroy and Li (2000).Mechanical turbulence controls atmospheric exchange during blowing snow, thus shear stress is calculated using the Prandtl-von Kármán logarithmic wind speed profile.
The aerodynamic roughness length z 0 (m) is controlled by the saltation height and is calculated using coefficients that were derived from extensive field observations in a level location with fully developed saltating flow (Pomeroy and Gray, 1990).Recent work has raised uncertainty in the turbulent exchange processes during blowing snow.Doorschot et al. (2004) showed that the aerodynamic roughness length in alpine terrain is more influenced by surrounding topography than by the saltation height.According to their findings, z 0 in this environment may be underestimated by the Pomeroy and Gray (1990) formulation.According to Helgason (2010) and Helgason and Pomeroy (2005), the advected turbulence associated with surrounding topography primarily enhances shear stress via horizontal turbulence and considerably less so via vertical turbulence.So the effects on particle lift and vertical fluxes such as sublimation may be muted.For this reason, it is assumed that snow-atmosphere exchange in alpine environments can be adequately simulated using roughness length calculations derived from measurements over saltation and vegetation roughness in open and flat terrain.
F susp is calculated as a vertical integration from a reference height near the top of the saltation layer, h * , to the top of the blowing snow boundary layer (z b ), given by Pomeroy and Male (1992) where k is von Kármán's constant (0.41), u * is the friction velocity, η is the mass concentration of blowing snow at height z (m) and z 0 is the aerodynamic roughness height.z b is governed by the time available for the vertical diffusion of snow particles from h * , calculated using turbulent diffusion theory and the logarithmic wind profile.h * increases with friction velocity and is estimated as given by Pomeroy and Male (1992).
For fully-developed flow it is constrained at z b =5 m.Shear stress is assumed constant (dτ /dt=0) up to z b and suspension occurs under steady-state conditions (dη/dt=0).As suspension diffuses from the saltation layer, saltation must be active for suspension to proceed.
E B is calculated as a vertical integration of the sublimation rate of a single ice particle with the mean particle mass being described by a two-parameter gamma distribution of particle size that varies with height.The calculation of the sublimation rate of a single ice sphere is based on Thorpe and Mason's (1966) relationships based on a balance of atmosphereparticle heat transfer and latent heat due to sublimation at the particle surface.The vertically integrated sublimation rate is given by where m is the mean mass of a single ice particle at height z.The rate that water vapour can be removed from the ice particle's surface layer, dm/dt, is calculated assuming particles to be in thermodynamic equilibrium.dm/dt is controlled by radiative energy exchange, convective heat transfer to the particle, turbulent transfer of water vapour from the particle to the atmosphere and latent heat associated with sublimation, and is given by Schmidt (1972) assuming the particle ventilation velocity to be equal to the particle terminal fall velocity.E B calculations are highly sensitive to ambient relative humidity, temperature and wind speed (Pomeroy et al., 1993;Pomeroy and Li, 2000).Field observations show that blowing snow is a phenomenon that is unsteady over both space and time.The time steps most frequently used in PBM studies (i.e. 15, 30 or 60 min) do not always match the highly variable and intermittent nature of blowing snow.In addition, small scale spatial variability in snowcover properties produces sub-element (e.g.grid cell or HRU) variability in snow transport.Li and Pomeroy (1997b) developed an algorithm to upscale blowing snow fluxes from point to area.The probability of blowing snow occurrence, p, is approximated by a cumulative normal distribution as a function of mean wind speed (location parameter), the standard deviation of wind speed (scale parameter).Empirical equations for the location and scale parameter were developed from six years of data collected at 15 locations in the Canadian prairies and are calculated from the number of hours since the last snowfall and the ambient atmospheric temperature.

SNOBAL
SNOBAL (Marks et al., 1998(Marks et al., , 1999) ) calculates the amount of snowmelt using the energy equation and was developed for deep mountain snowcover.SNOBAL represents the snowcover as two layers: a surface active layer of fixed depth and a lower layer that represents the remaining snowpack.SNOBAL was previously applied at MCRB by DeBeer and Pomeroy (2009).
Snowmelt in either layer occurs when the available energy exceeds that required to bring the snow layer temperature above 0 • C. The amount of snowmelt (depth of water equivalent per unit time) is calculated using where Q m is the energy available for melt, ρ is the density of water, h f is the latent heat of fusion (333.5 kJ kg −1 ) and B is the fraction of ice in a unit mass of wet snow (0.97).Q m is calculated from the energy equation as where Q n is the net radiation (incoming and outgoing shortwave and longwave radiation), Q h is the convective sensible heat flux, Q e is the convective latent heat flux, Q g is the conductive ground heat flux, Q p is the advective heat from rainfall, Q A is the small-scale advective heat from bare ground and U/ t is the change in internal energy of the snow mass.Snow albedo, α, during winter is estimated using the method outlined by Gray and Landine (1987).The albedo depletion was approximated by three lines of different slope representing three periods: pre-melt, melt and post-melt.
Sublimation (and condensation), E, at the snowatmosphere interface is diagnosed from the latent heat flux, and sublimation (condensation) at the snow-soil interface is calculated from the specific humidities of the snow and soil layers and a diffusion coefficient (Marks et al., 1998).Liquid water in the snowpack is first evaporated using the ratio of latent heat of vaporization to sublimation at 0 • C (0.882).The remaining diagnosed evaporative loss is calculated as ice sublimation.Half of the sublimated ice decreases snowpack depth but does not alter its density.All of the evaporated liquid water and the other half of sublimated ice decrease the snowpack density and specific mass.
SNOBAL ensures numerical stability by using a variable time-stepping scheme that is dependent on the snowpack specific mass.There are three snow mass thresholds (60, 10 and 1 kg m −2 or mm SWE) that activate the use of different time step durations.The longest time step is used when the snow mass is above 60 kg m −2 and the shortest time step is used when the snow mass is below 1 kg m −2 .

Canopy
The CRHM Canopy module calculates energy and water input to the snow surface beneath a needleleaf forest canopy and is fully described in Ellis et al. (2010).
Shortwave transmissivity through the forest canopy is given by a Beer-Bouger type formulation given by a variation of Pomeroy and Dion's (1996) formulation (Pomeroy et al., 2009) where θ is the solar angle above the horizon (radians) and PAI' is the effective winter plant area index.Multiple reflections within the canopy are not explicitly account for and there are not separate calculations for canopy foliage, trunks and gaps.Enhanced longwave irradiance to the surface from the forest canopy is calculated.The incoming longwave radiation to the snowpack beneath the forest canopy, L ↓f , is given by where v is the sky view factor, ε f is the forest thermal emissivity (unitless), σ is the Stefan-Boltzmann constant (W m −2 K −4 ) and T f is the forest temperature (K).
Canopy also estimates the canopy throughfall of rain and snow, the canopy interception and evaporation of rain, the canopy interception and sublimation of snow, the unloading of intercepted snow and the drip of intercepted rain.
Canopy interception of snowfall and sublimation of intercepted snow is calculated using relationships presented by Pomeroy et al. (1998).The amount of intercepted snow is calculated as where C 1 is the dimensionless canopy-leaf contact per ground, P S is snowfall and I * s is the maximum intercepted snowload which is estimated as a function of the maximum snowload per unit area of branch, the density of falling snow and effective leaf area index (LAI').The sublimation of intercepted snow is calculated by adjusting the sublimation rate of an ice sphere by an intercepted snow exposure coefficient (Pomeroy and Schmidt, 1993).

Fisera ridge parameterization
HRUs were selected by grouping snow depths measured along the FR transect (Fig. 4).These manual snow depth measurements captured the spatial variability in wind exposure along the FR transect, which exerts a stronger control on winter snow accumulation at this location than does solar radiation and vegetation.Observed snow depth from the five dates shown on Fig. 4 show that the spatial snow depth pattern is temporally stable over winter.Five HRUs were selected based on the observed snow depth patterns shown in Fig. 4. The north-facing slope HRU (NF) is located from 127 to 181 m, the ridgetop HRU (RT) is located from 90 to 127 m, the upper southfacing slope (SF-upper) is located from 28 to 90 m, the lower south-facing slope (SF-lower) is located from 0 to 28 m and the Forest HRU is located from 0 to −15 m.
The HRUs follow a fixed aerodynamic sequence in that the model always transports snow from upwind to downwind HRUs.The HRU snow transport sequence is NF→RT →SFupper → SF-lower→Forest (i.e.NF snow transport reaches all of RT, SF-upper, SF-lower and Forest; SF-upper snow transport only reaches SF-lower and Forest; etc.).The static definition of the HRU locations and relative lengths is a simplified representation of the actual spatiotemporal snow redistribution patterns.This assumption is appropriate at this location.Snow transport rates scale approximately with the fourth power of wind speed (Pomeroy and Male, 1992;Essery et al., 1999).Fig. 5 shows a fourth power of wind speed (u 4 ) rose (sum of u 4 binned by direction) for the observed wind speed and direction at the Fisera Ridge ridgetop station over the study periods.It is clear from this figure that northerly-to-northwesterly blowing snow events dominate snow redistribution at this location.Furthermore, changing HRUs sizes during a model run would add considerable complexity to the calculation of mass balances for HRUs.
Key CRHM model parameters are presented in Table 2.The Canopy module was only applied to the Forest HRU.All CRHM model parameters were set based on either field measurements or default/typical values with the exception of vegetation height on the NF and RT.Shorter shrub heights than measured were needed to scour enough snow from these HRUs.PBSM is parameterized for densely spaced, narrow crop stalks and grass.Shorter vegetation heights parameters were required to represent sparse shrubs on the NF and RT HRUs.This indicates that the PBSM parameterization for the aerodynamic roughness height of vegetation may not be suitable for such shrubs and should be revised.Average HRU aspect and slope were determined from the DEM.A blowing snow fetch distance of 300 m was specified for each HRU as this is the minimum value required for the fully-developed flow calculations performed by PBSM.For SNOBAL calculations, the maximum active layer thickness was fixed at 0.1 m and snow surface roughness length was set to 0.001 m.Simulations were performed for 2007/2008 and 2008/2009 applying the ridge-top station air temperature, relative humidity and incoming longwave radiation observations to all HRUs.The north-facing meteorological station wind speed data was applied to the NF, the ridge-top meteorological station wind speed data was applied to the RT, and the southeast-facing meteorological station wind speed was applied to the SF-upper, SF-lower and Forest.Incoming shortwave radiation observations from the ridge-top station (considered a flat plane) were applied to each HRU after adjustments for aspect and slope made by the Global and Slope Qsi modules.Reflected radiation from adjacent terrain was captured by the radiometer measurements at the ridge-top station, therefore all HRUs received the same contribution of reflected radiation relative to total incoming radiation.This approach for taking into account reflections deviates somewhat from reality; however this approach produced excellent radiation for snowmelt modelling in the same environment (DeBeer and Pomeroy, 2009).It is not necessary to model topographic shading for this study because it is accounted for in the ridge-top measurements and it is suitable to assume identical effects over this short model transect.
Two simulations were performed for each period using the same parameterization in order to investigate the suitability of PBSM's sublimation calculations in this environment and the importance of blowing snow sublimation on the snow mass balance in this environment: 1. a baseline simulation that includes calculations of blowing snow transport and sublimation, snowpack sublimation and snowmelt; and 2. a simulation omitting the calculation of blowing snow sublimation.
The suppression of blowing snow sublimation calculations has a non-linear effect on snow transport rate calculations.
Given a surface snow erosion rate, it is obvious that proportionally more snow is transported when sublimation is suppressed.However, the steady-state snow erosion rate is partly limited by the instantaneous drift density.Thus, the surface snow erosion rate is concomitantly suppressed with the blowing snow sublimation rate.3 and 4 show end-of-winter snow accumulation, cumulative snowmelt, transport in to and out of HRUs, blowing snow sublimation, snowpack sublimation and sublimation of intercepted snow in the Forest HRU for the baseline simulations and simulations omitting blowing snow sublimation, respectively.underestimated throughout the simulation, respectively, and describes the reproduction of total snow mass over all the HRUs.The RMSE gives a measure of the variation of residuals between observed and simulated snowcover in mm SWE and described how well the snow mass is distributed on the various HRUs.The coefficient of determination gives a measure of the accuracy of the model in simulating the observed SWE.A R 2 =1.0 indicates that the model perfectly simulated the variation in observed SWE and R 2 =0.0 indicates that the model did not simulate any of the variation.Snow accumulation was well simulated with CRHM when using the observed wind speed data.2007/2008 was not simulated to the same accuracy as 2008/2009 as snow accumulation was overestimated on the SF-upper throughout the simulation and snow accumulation was overestimated on the SF-bottom during February but matched observed snow Omitting blowing snow sublimation calculations resulted in an overestimation of snow accumulation on all HRUs during both periods.Though more snow transport was simulated when blowing snow sublimation calculations were omitted, the total amount of snow eroded from the surface was lower because it was limited by the ability of the atmosphere to transport snow.Omitting blowing snow sublimation calculations resulted in a substantial overestimation of the snow mass available for melt at the end of winter (Table 5).This error would cause substantial difficulties in accurately simulating snowcover ablation and runoff during melt-dominated periods from May-June.

Discussion
Estimated blowing snow sublimation losses over the transect were 19% for 2007/2008 and 17% for 2008/2009.These blowing snow sublimation losses were substantial and important to the winter water balance of the alpine ridge.Sublimation losses of this magnitude are considered to be a conservative estimate of alpine blowing snow losses in the Canadian Rocky Mountains since the model transect is located in a low alpine area where windflow separation over the ridge-top was not observed.It is believed that blowing snow sublimation losses are greater in the high alpine zone where the topography is sharper and steeper and flow separation occurs.
Satisfactory FR snow mass balance closure suggests that PBSM can be applied with caution in this environment and that the use of the minimum PBSM fetch distance parameter (300 m) is adequate in this environment.Boundary layer development for fetches shorter than this in complex terrain are poorly understood and so the parameter is left to its minimum value (based on the limits of PBSM physics) until a more realistic parameterization is be developed.
The observed SF-lower snow accumulation was greater than the Forest snow accumulation in 2008/2009, whereas the opposite was true during 2007/2008 and for the simulations.Observed wind speeds were generally higher during 2007/2008 than 2008/2009 (higher mean and less positive skew of wind speed), so it may have been a case of downwind transport distance increasing with increasing wind speed.

Conclusions
Snow redistribution and sublimation by wind, snowpack sublimation, snowmelt and the resulting accumulation regimes were simulated over HRUs representing a transect along an alpine ridge in the Canadian Rockies.This study shows that snow accumulation can be adequately simulated in HRUs over mountainous terrain using physically based blowing snow and snowcover models.A HRU-based discretization can be a much more computationally efficient approach than a fully-distributed one.This is particularly relevant for modelling snow redistribution within large-scale hydrology models and land surface schemes.HRUs were selected by examining manual snow depth measurements.Future work will involve generalizing HRUs based on terrain characteristics.
Snow mass losses due to blowing snow sublimation and snowpack sublimation were shown to be significant in this environment.Snow transport and sublimation from windward slopes and ridge-tops reduced snow accumulation in these landscapes to approximately one-quarter of snowfall.Approximately half of this loss can be attributed to blowing snow sublimation and snowpack sublimation, of which three-quarters was due to blowing snow sublimation alone.Transect blowing snow sublimation and snowpack sublimation losses over the winter periods were estimated to be 17-19% and 1-15% of cumulative snowfall, respectively.The estimated snow mass loss due to blowing snow sublimation is considered to be a conservative estimate for the alpine zone of the Canadian Rocky Mountains since windflow separation was not observed over this ridge.
An investigation into the importance of estimating blowing snow sublimation on the simulated snow mass balance was made by performing identical model runs but neglecting the calculation of blowing snow sublimation.Snow accumulation was considerably overestimated throughout the winter and end-of-winter snow accumulation was overestimated by approximately 30%.This would cause further difficulties in accurately simulating snowcover ablation and runoff during snowmelt-dominated periods.These results show that a physically based blowing snow model developed from measurements made in open flat areas can be applied over HRUs in mountainous terrain.
developed a snow drift index to predict lee slope drifting over four terrain aspects.MacDonald et al. (2009) built upon Dornes et al.'s work and identified hydrological response units (HRUs) suitable for calculating snow redistribution calculations in subarctic mountains with moderate topographic roughness and ran a blowing snow model to estimate snow accumulation quantities that compared well to field measurements.

Fig. 3 .
Fig. 3.Control volume for blowing snow mass fluxes and snowmelt energy.

Fig. 4 .
Fig. 4. Observed snow depth along Fisera Ridge transect.The boundaries of the five HRUs are indicated by the dashed lines.

Table 3 .
Summary of cumulative baseline model results for (a) 2007/2008 and (b) 2008/2009 (quantities are in kg m −2 ; brackets indicate quantity as percentage of snowfall).

Table 4 .
Summary of cumulative model results without blowing snow sublimation for (a) 2007/2008 and (b) 2008/2009 (quantities are in kg m −2 ; brackets indicate quantity as percentage of snowfall).

Table 5 .
Model evaluation statistics for baseline simulations and simulations omitting blowing snow sublimation.

Table 6 .
Model evaluation statistics only for final pre-melt measurement dates for baseline simulations and simulations omitting blowing snow sublimation.