1430

JOURNAL OF HYDROMETEOROLOGY

VOLUME 10

A Physically Based Daily Hydrometeorological Model for Complex Mountain Terrain RYAN J. MACDONALD, JAMES M. BYRNE, AND STEFAN W. KIENZLE Department of Geography, University of Lethbridge, Lethbridge, Alberta, Canada (Manuscript received 10 September 2008, in final form 29 June 2009) ABSTRACT This paper describes the continued development of the physically based hydrometeorological model Generate Earth Systems Science input (GENESYS) and its application in simulating snowpack in the St. Mary (STM) River watershed, Montana. GENESYS is designed to operate a high spatial and temporal resolution over complex mountainous terrain. The intent of this paper is to assess the performance of the model in simulating daily snowpack and the spatial extent of snow cover over the St. Mary River watershed. A new precipitation estimation method that uses snowpack telemetry (SNOTEL) and snow survey data is presented and compared with two other methods, including Parameter-elevation Regressions on Independent Slopes Model (PRISM) precipitation data. A method for determining daily temperature lapse rates from NCEP reanalysis data is also presented and the effect of temperature lapse rate on snowpack simulations is determined. Simulated daily snowpack values compare well with observed values at the Many Glacier SNOTEL site, with varying degrees of accuracy, dependent on the method used to estimate precipitation. The spatial snow cover extent compares well with Moderate Resolution Imaging Spectroradiometer (MODIS) snow cover products for three dates selected to represent snow accumulation and ablation periods.

1. Introduction Water supply in western North America is dependent on snowpack from mountainous regions (Barnett et al. 2005; Field et al. 2007; Mote et al. 2005). The complex interaction between snowpack and meteorological variability makes these regions extremely vulnerable to changes in climatic processes (Beniston 2003; Leung and Wigmosta 1999; McKenzie et al. 2003). It is expected that mountain snow accumulations will decline with continued atmospheric warming (Hamlet and Lettenmaier 1999), resulting in a reduction of available water from snowpack (Barnett et al. 2005; Lapp et al. 2005). Mountain snowpack plays an important role in almost every component of the hydrological balance. Snow cover has an effect on local meteorological conditions, soil moisture conditions (Groisman et al. 1994; Kane et al. 1991; Zhang et al. 2003), the distribution and growth season of vegetation (Stephenson 1990), and the timing and availability of runoff (R; Fontaine et al. 2002). The importance of snow in mountainous regions has led to significant Corresponding author address: James M. Byrne, Department of Geography, University of Lethbridge, 4401 University Drive, Lethbridge, AB T1K 3M4, Canada. E-mail: [email protected] DOI: 10.1175/2009JHM1093.1 Ó 2009 American Meteorological Society

research in snow hydrology and the development of spatial hydrometeorological models. Hydrometeorological measurements in mountainous regions are sparse (Marks et al. 1992) and do not represent the variability required for modeling entire watersheds (Diaz 2005). Spatial estimates of hydrometeorology in mountainous environments are, therefore, frequently made using a low number of point measurements as input to spatial models (Liston and Elder 2006b). Spatial models rely on the interaction between physiographic characteristics of the landscape and meteorological processes to make estimates of hydrometeorological variables. In mountainous environments, the high spatial and temporal variability in hydrometeorological conditions requires spatial models that are physically realistic and computationally efficient (Liston and Elder 2006b). A number of models have been developed to simulate mountain hydrometeorology. The Regional Hydroecological Simulation System (RHESSys; Band et al. 1991, 1993), the Precipitation–Runoff Modeling System (PRMS; Leavesley et al. 1983), snow evolution model (SnowModel; Liston and Elder 2006a), and Alpine3D (Lehning et al. 2006) are four models that have been used for hydrological and ecological modeling in mountainous watersheds. RHESSys integrates GIS and a series of

DECEMBER 2009

MACDONALD ET AL.

1431

FIG. 1. The STM River watershed in Montana and southern Alberta.

subprograms to spatially estimate ecosystem processes at the watershed scale. PRMS is a distributed model that was designed to evaluate the effects of precipitation ( p), climate, and land use on general basin hydrology, whereas SnowModel is a detailed spatial snowpack model developed for application under a range of landscapes where snow occurs. Alpine3D is a surface energy balance model that has been used to simulate finescale snow processes in mountainous regions. Although these models are useful, they are not always available and options may be limited for application to large watersheds with limited data. This paper describes the continued development and application of a model for simulating hydrometeorological conditions in large watersheds with relatively little observed data. Generate Earth Systems Science (GENESYS) is a physically based model for spatially estimating daily hydrometeorological variables over mountainous terrain using routinely available meteorological data. GENESYS is under development at the University of Lethbridge under the direction of Dr. James Byrne and has been applied in several studies (A. Sheppard 1996, personal communication; Lapp et al. 2002, 2005; Larson et al. 2009, manuscript submitted to J. Hydrol., hereafter LBJLK), including work described herein. The objective of this paper is to further develop the GENESYS model to more accurately represent mountain hydrometeorology spatially and determine how well the model simulates snow accumulation and ablation.

2. Study area and meteorological data The headwaters of the St. Mary River watershed lie on the eastern slopes of the Rocky Mountains, with the majority of the upper watershed residing within Glacier National Park, Montana. The St. Mary River flows from the continental divide, through the upper and lower St. Mary lakes, and ends in southern Alberta, where it meets the Oldman River (Fig. 1). The climatic regime is a transitional zone between coastal and continental climates. The region is also influenced by the orographic effect, which is most noticeable during the winter months because synoptic conditions dominate (Hanson 1982). The area receives the majority of its precipitation in the winter, with snow accounting for roughly 70% of the annual precipitation at high elevations (Selkowitz et al. 2002). The total drainage area of the study watershed is 1195 km2, with a mean elevation of 1745 m, ranging from 1249 to 3031 m. The area is a relatively undisturbed, ecologically diverse region, which is largely attributed to the fact that a large portion of the drainage area is within Glacier National Park. Coniferous forests account for 24% of the land cover, deciduous forests account for 21%, and herbaceous plants cover another 29% of the area; 23% of the area is barren rock or soil, and 3% of the area is water (USGS 2006). The St. Mary climate station was selected to drive the model. This station is centrally located at an elevation of 1390 m near St. Mary, Montana, in the eastern portion

1432

JOURNAL OF HYDROMETEOROLOGY

of the watershed. Daily temperature and precipitation data for the period from 1960 to 2005 at the St. Mary climate station were obtained (NCDC 2006). There were significant data gaps in the station record from 1960 to 1982, with minor data gaps from 1982 to 2005. Therefore, missing records from the years 1960 to 2005 were infilled using nearby climate stations and linear regression (LBJLK). To derive precipitation–elevation relationships, snow water equivalent (SWE) measurements from the Preston snow survey (Fig. 1) were used. The snow survey is operated by the U.S. Geological Survey (USGS); it began in 1994 and continues to the present. The survey has 32 sampling points located near the center of the watershed and spans an elevation range from 1438 to 2290 m. SWE data have been acquired from the inception of the survey to the end of the 2006 snow year (D. Fagre 2006, personal communication). SWE data from the Many Glacier (MG) snowpack telemetry (SNOTEL) site (Fig. 1) have also been used in this study. The site is located in a small meadow surrounded by trees at an elevation of 1519 m in the western portion of the basin (NRCS 2006). This site has been in operation since 1976 and continues to the present; daily SWE data were obtained for the period from 1976 to 2005.

3. Model description The GENESYS model is designed to operate at high spatial and temporal resolution using two individual components to simulate the hydrological balance over mountainous terrain. The first component, SimGrid (A. Sheppard 1996, personal communication), applies the Mountain Climate Simulator (MTCLIM) model (Hungerford et al. 1989) and GIS-derived modeling units referred to as terrain categories (TCs). The MTCLIM model is looped in SimGrid to provide daily estimates of temperature, precipitation, solar radiation, and relative humidity for the subsequent modeling of hydrological processes in each TC carried out by the SnowPack component of the model. The SnowPack component applies physical equations to simulate sublimation, canopy interception, snowmelt, soil water storage, evapotranspiration, and runoff from each TC.

a. Derivation of terrain categories The complex nature of mountainous environments implies that high-resolution spatial data are needed for hydrometeorological simulations to be relevant. However, it is important to understand that increased resolution implies increased complexity and not necessarily a higher degree of accuracy (Daly 2006). Therefore, the explana-

VOLUME 10

tion of hydrometeorological variables at an appropriate spatial scale in mountainous terrain is important. Given the physical structure of the MTCLIM model, we suggest that an appropriate spatial scale can be determined using ecosystem responses to climatic processes. Because vegetative cover is highly dependent on hydrometeorological conditions (Mather and Yoshioka 1968; Stephenson 1990), it provides an ecologically sensitive surrogate for the spatial variability in hydrometeorological conditions. A combination of vegetative cover and elevation is used to represent spatial variability in hydrometeorology over the St. Mary River watershed. A land cover grid derived using Landsat imagery (USGS 2006) was overlaid with a 100-m digital elevation model (DEM) classified into 100-m elevation intervals to determine TCs for the St. Mary basin. The land cover grid consisted of nine categories: dry herbaceous, mesic herbaceous, deciduous trees–shrubs, coniferous trees– shrubs, coniferous trees–open, water, snow, barren rock–soil, and shadows. Snow and shadow classes were eliminated from the land cover grid and assigned the values of the nearest land cover. The combination of elevation and land cover resulted in 82 TCs over the St. Mary River watershed. TCs range in area from 100 m2 in the topographically heterogeneous portions of the watershed to 88 km2 in the low elevation, relatively homogenous portions of the watershed. For each TC mean slope, aspect, and elevation values were derived. The MTCLIM model was applied to all 82 TCs.

b. Application of MTCLIM The MTCLIM model uses two types of climatological logic: a topographic logic that determines meteorological variables by extrapolating data from a base climate station to the TC, and a diurnal climatology that derives additional information from climate station data (Hungerford et al. 1989). The diurnal climatology in MTCLIM generates incident solar radiation and relative humidity, whereas the topographic logic extrapolates climate station data to make estimations of maximum and minimum air temperature and precipitation (Glassy and Running 1994). MTCLIM can be driven by any climate station that provides maximum and minimum temperatures and precipitation. Climate station data used to drive the MTCLIM model can be referred to as base data. For each TC, MTCLIM requires mean elevation, mean slope, mean aspect, mean monthly precipitation, and monthly-mean leaf area index (LAI) values from the 1-km Moderate Resolution Imaging Spectroradiometer (MODIS)/Terra global dataset (Roy et al. 2002). Variables set as constants over all TCs are surface albedo (0.2), and atmospheric transmissivity (0.65). Here only the precipitation and temperature methods

DECEMBER 2009

MACDONALD ET AL.

1433

FIG. 2. Seasonal distribution of precipitation at MG and STM.

within MTCLIM are reported. For a more detailed description of relative humidity and solar radiation estimates, refer to Hungerford et al. (1989) and Glassy and Running (1994).

1) PRECIPITATION Of the hydrological variables, precipitation is perhaps the most difficult to quantify spatially. By accounting for the physiographical controls on spatial and temporal distribution of precipitation and incorporating observed meteorological data, estimates of precipitation can be made for mountainous environments (e.g., Daly et al. 2008). However, because of the complexities in processes controlling precipitation in mountainous terrain, it is difficult to estimate precipitation at the daily time step. At a coarser temporal scale, the variability in precipitation is reduced and, therefore, it can be better described. Applying monthly data to adjust daily precipitation values enables large-scale precipitation patterns to be maintained while accounting for daily variability. A similar method is used by Running et al. (1987), where they apply annual data. However, determining monthly precipitation values over the entire watershed is difficult. An objective of this study is to determine the most suitable monthly spatial precipitation estimation method by comparing how precipitation inputs affect daily SWE

simulations at the Many Glacier SNOTEL site (NRCS 2006). Three methods are applied; two of the methods use observed data within the St. Mary River watershed to derive precipitation–elevation relationships. The first method was developed by LBJLK, where precipitation estimates are made using a precipitation–elevation function. The second method is presented here, where precipitation estimates are made as a function of elevation and season. The third method applies the 1971–2000 monthly precipitation averages from the Parameterelevation Regressions on Independent Slopes Model (PRISM) dataset (Daly et al. 2008). The climatological characteristics of the St. Mary and Many Glacier sites were assessed. Figure 2 shows that the seasonal distribution of precipitation differs significantly between the St. Mary and Many Glacier sites (r 2 5 0.09, p 5 0.17). Three methods are presented to show the relevance of accounting for differences in seasonality between low-elevation climate stations and mountainous regions of the watershed. One method is shown that does not account for this seasonal difference, whereas the other two methods are shown that do account for seasonality.

(i) Precipitation method A LBJLK derived a method that established precipitation– elevation relationships. The method applies a linear

1434

JOURNAL OF HYDROMETEOROLOGY

VOLUME 10

FIG. 3. Linear relationship between local elevation in relation to the MG SNOTEL site and dSWE at the Preston snow course.

precipitation–elevation function to daily data [Eq. (1)]. The equation was derived using monthly changes in SWE at the Preston snow course: P(Daily) 5 Pstm(Daily) 1 0.232 3 elevation 3

pstm(Daily) pstm(Monthly)

,

(1) where P(Daily) is daily precipitation (mm) at the TC, Pstm(Daily) is daily total precipitation (mm) at St. Mary, Pstm(Monthly) is monthly precipitation averages (mm) for St. Mary, and elevation (m) is the local elevation relative to St. Mary, where the elevation at St. Mary is set to 0 m. This method is applied in the first run of the GENESYS model and assessed using daily SWE data from the Many Glacier SNOTEL site.

(ii) Precipitation method B The second method applies MTCLIM logic (Running et al. 1987), using monthly data to adjust daily precipitation values as a function of elevation over the watershed. This is done by calculating a ratio between mean monthly precipitation values at each TC and the St. Mary climate station. The monthly ratios are multiplied by the daily precipitation value to adjust daily data over the watershed. Data used to derive this method included the Preston snow survey, the St. Mary climate station, and the Many Glacier SNOTEL site. The following series of steps was used to derive precipitation values at each TC: 1. Derive a precipitation–elevation relationship. 2. Adjust for seasonality differences between the St. Mary climate station and the mountain portions of the watershed. 3. Calculate ratios between monthly-mean values at the St. Mary climate station and each TC. 4. Apply ratios to daily precipitation data.

A change in winter SWE (DSWE) was calculated for each monthly sampling interval for 73 months (LBJLK). These monthly values were only calculated for the cold season (from January to March), where average temperatures for the sampling period were below 08C. Negative DSWE values were omitted from the analysis, with the assumption that these values corresponded to melt. With these two constraints, it is assumed that monthly DSWE values correspond with monthly increases in precipitation. The two precipitation–elevation relationships are derived relative to the St. Mary climate station and the Many Glacier SNOTEL site. Using precipitation–elevation relationships relative to both St. Mary and Many Glacier enables the model to account for seasonality differences between mountain and low land areas. Local elevations are determined relative to both St. Mary and Many Glacier, where St. Mary and Many Glacier were set to have a local elevation of 0 m. The resultant relationships between local elevation and mean winter DSWE are presented in Figs. 3 and 4 . An adjustment was made to account for the seasonal differences in precipitation between the St. Mary climate station and the Many Glacier SNOTEL site, resulting in a better representation of precipitation over the watershed. Monthly relationships between St. Mary and Many Glacier precipitation means for the years 1982–2005 were derived using linear regression (Table 1). The time period from 1982 to 2005 was selected because of data gaps prior to 1982 at the St. Mary climate station. The St. Mary climate station best represents the Many Glacier SNOTEL site during the winter months. On the basis of the regression results listed in Table 1, it is assumed that the St. Mary climate station can be reliably used to predict monthly precipitation at Many Glacier during the winter. The transitional months of April, May, and September had relatively poor monthly precipitation relationships. However, because of the limited available

DECEMBER 2009

1435

MACDONALD ET AL.

FIG. 4. Same as Fig. 3, but for STM and dSWE.

data, all 12 monthly relationships were applied to predict a mountain base station, which in this case is the Many Glacier SNOTEL site. This was done to enable the model to account for changes in seasonality using St. Mary as the single low-elevation base station, where the longest climate records are available. To determine monthly precipitation means at each TC, the two DSWE equations are used. For TC elevations below 1500 m, the Preston DSWE relationship relative to St. Mary (Fig. 4) is applied with the monthlymean values from St. Mary climate station as the base: SWE(,1500) 5 0.051(local elevation) 1 base.

(2)

If the TC elevation is 1500 m and above, monthlymean precipitation values for a mountain base are calculated using monthly relationships presented in Table 1. Using the predicted mountain base, the Preston DSWE relationship relative to Many Glacier (Fig. 3) is applied: SWE(.1500) 5 0.063(local elevation) 1 mountain base. (3) This allows for a seasonal shift in precipitation to be made for the mountainous portion of the watershed, while accounting for the effects of elevation. At elevations greater than 2300 m (extent of the snow course data), monthly means are assigned the same value as the mean at 2300 m, resulting in no change in SWE with elevation above 2300 m. Figure 5 shows how simulated monthly-mean precipitation volumes change with season and elevation. To apply the effect of elevation and shift in seasonality to the daily historical precipitation record, ratios are calculated between monthly precipitation means at the St. Mary climate station and the monthly means at each TC. This results in the largest ratios during the winter and smallest ratios during the summer (Fig. 6). These ratios are multiplied by daily precipitation volumes at the St. Mary climate station, resulting in pre-

cipitation volumes that are adjusted as a function of elevation and season over the watershed.

(iii) Precipitation method C The third precipitation estimation method applies the 1971–2000 precipitation means from PRISM to obtain monthly precipitation values. For each of the 12 monthly surfaces, a mean precipitation value is calculated for each TC. To maintain consistency at the 1400-m elevation band, precipitation values are used from the St. Mary climate station. The monthly values derived are used to determine the monthly ratios between the St. Mary climate station and each TC within the watershed. Figure 7 demonstrates that PRISM accounts for the differences in seasonality between St. Mary and the higher mountainous portions of the watershed. It is important to note that PRISM is also able to account for other topographic influences on precipitation, such as coastal proximity and slope orientation (Daly et al. 2002). To maintain consistency between comparisons, the TCs shown in Fig. 7 are the same as the TCs shown in Fig. 6. Using PRISM inputs results in slightly higher ratios when compared to method B. However, they still reflect TABLE 1. Linear relationships between mean monthly precipitation at STM and mean monthly precipitation at MG (n 5 23). Month

Mountain base equation

r2

p

Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec

MG 5 1.526(STM) 1 58.575 MG 5 1.505(STM) 1 29.284 MG 5 1.587(STM) 1 28.825 MG 5 1.233(STM) 1 28.189 MG 5 0.549(STM) 1 55.335 MG 5 0.863(STM) 1 33.151 MG 5 0.735(STM) 1 26.554 MG 5 0.875(STM) 1 13.842 MG 5 0.904(STM) 1 31.279 MG 5 1.875(STM) 1 17.199 MG 5 1.797(STM) 1 45.295 MG 5 1.570(STM) 1 46.428

0.83 0.77 0.82 0.57 0.29 0.83 0.85 0.72 0.53 0.84 0.89 0.87

,0.0001 ,0.0001 ,0.0001 ,0.0001 0.006 ,0.0001 ,0.0001 ,0.0001 ,0.0001 ,0.0001 ,0.0001 ,0.0001

1436

JOURNAL OF HYDROMETEOROLOGY

VOLUME 10

FIG. 5. Mean monthly precipitation change as a function of season and elevation.

the change in seasonality between mountainous portions of the watershed and the St. Mary climate station. These ratios are used to adjust the daily precipitation volumes at the St. Mary climate station as a function of elevation and season over the watershed.

2) TEMPERATURE To account for temperature changes as a function of elevation, MTCLIM applies temperature lapse rates. This study involves varying temperature lapse rates for three separate model runs to determine the effect of lapse rate on simulated snow accumulation and ablation. The first model run applies lapse rates of 8.28C km21 for maximum temperature and 3.88C km21 for minimum temperature. These lapse rates were used by LBJLK and resulted in temperature estimates that compared very well to an alpine site on Lakeview Ridge near Waterton, Alberta. The second run applies lapse rates derived at Castle Mountain ski resort, approximately 100 km northwest of the St. Mary River watershed by Pigeon and Jiskoot (2008). They determined maximum and minimum

temperature lapse rates to be 6.18 and 5.98C km21, respectively. The third model run applies daily lapse rates from 1961 to 2000 derived from National Centers for Environmental Prediction (NCEP) reanalysis data (NCEP 2008). To derive daily lapse rates, a method is used that calculates the differences in elevation and temperature between the 1000- and 700-mb surfaces and from those differences derives linear lapse rates (D. Blair 2008, personal communication). The four grid cells covering the watershed were selected and averaged for this analysis. The 1961–2000 average of daily NCEP lapse rates were 6.58 and 4.68C km21 with a range of 17.88 and 17.68C for maximum and minimum temperatures, respectively, accounting for inverted lapse rates.

c. Snowpack Daily spatial hydrometeorological data are used to model changes in TC snow water equivalent (SWE). When snowpack is present, a daily hydrological balance is calculated according to Eq. (4):

FIG. 6. Varying monthly-mean precipitation ratios between STM and TCs at the elevation bands of 1400, 2000, and 2800 m.

DECEMBER 2009

1437

MACDONALD ET AL.

FIG. 7. Same as Fig. 6, but using PRISM data.

SWE(t) 5 SWE(t1) 1 P(t)  I (t)  S(t)  R(t)  IF(t) , (4) where, SWE is the amount of snow water equivalent (mm) in the snowpack, P is simulated daily precipitation as rain or snow, I is canopy interception, S is sublimation, IF is infiltration, and t is the time step (days). If the snowpack has completely melted, a hydrological balance is calculated that accounts for evapotranspiration (ET) and changes in soil moisture (SM) conditions: SM(t) 5 SM(t1) 1 P(t)  I (t)  ET(t)  R(t) ,

(5)

which allows for two separate hydrological balances to be calculated.

1) PRECIPITATION PARTITIONING

(i) Snow interception

The method used to partition rain and snow was developed by Kienzle (2008). This method uses an S-shaped curve and two temperature variables:     T  Tt 2 T  Tt 2 Pr 5 5 1 6.76 1.4Tr 1.4Tr   T  Tt 1 3.19 1 0.5, 1.4Tr

estimates of canopy load, interception, and unloading. This model is adopted as a subroutine in the GENESYS model. This adaptation assumes that the cold climatic regime and physical properties of tree species of the boreal forest represent the east slopes of the Rocky Mountains well enough to provide realistic estimates of snow interception. An empirical rainfall interception routine based on canopy LAI was adapted for GENESYS (Von Hoyningen-Huene 1983). Snow interception (snowInt; mm SWE) and rain interception (rainInt; mm) by the canopy are determined for each terrain category containing coniferous forests. Only rain interception is calculated for deciduous forests because a subroutine was not available that accounts for snow interception in deciduous forests.

The snowInt subroutine uses the following formula derived by Hedstrom and Pomeroy (1998): SnowInt 5 I 3 0.678,

(6)

where Pr is the proportion of precipitation that falls as rain (range from 0 to 1), T is the mean daily temperature, Tt is the threshold mean daily temperature, and Tr is the range of temperatures where both rain and snow can occur. Values for Tt and Tr of 3.78 and 168C, respectively, are used, as suggested by Kienzle (2008) for a similar climate station in southern Alberta.

2) CANOPY INTERCEPTION The canopy snow interception model developed by Hedstrom and Pomeroy (1998) for the southern boreal forests of western Canada provides physically based

(7)

where 0.678 is a determined unloading coefficient for pine and spruce forests and I is the intercepted snow load at the start of unloading (Hedstrom and Pomeroy 1998); I is determined using the following formula: I 5 (L  Load)(1  10kP ),

(8)

where Load is the total snow load in the canopy on the previous day (kg m22), L is the maximum load for the given canopy given the boundary layer conditions (kg m22), and P is the amount of precipitation falling as snow. An approach is taken for the determination of k (the proportionality factor), in which it is assumed that there is a closed canopy and interception is completely efficient. It is also assumed that snowflakes are falling vertically on the canopy. The following formula represents k when the preceding assumptions are made:

1438

JOURNAL OF HYDROMETEOROLOGY

k5

1 , L

(9)

where L 5 S 3 LAI.

(10)

LAI values were obtained for each TC from MODIS, and S is the species value at a given density; S is determined using the following formula:    46 . S 5 SV 0.27 1 Ps

(11)

The value SV is a constant equal to 5.9 kg m22 derived by Hedstrom and Pomeroy (1998) for spruce forests; Ps (kg m23) is the density of snow calculated as a function of mean air temperature (Hedstrom and Pomeroy 1998). The canopy load is calculated by adding each interception event to the canopy store and subtracting snow that is sublimated. The canopy is able to store snow until L is reached, at which time the remaining snow in the canopy will fall to the ground and is then incorporated in the snowpack.

(ii) Rain interception RainInt is estimated using the Von Hoyningen-Huene (1983) formula that calculates interception as a function of total rainfall and LAI. RainInt is calculated on a daily basis using the following function: RainInt 5 0.30 1 0.27Rain 1 0.13LAI  0.013Rain2

in forested areas to be minimal. The sublimation model requires three environmental inputs calculated by the SimGrid subroutine: total incident radiation Q (W m22), * daily average temperature Ta (K), and relative humidity RH (%). Total sublimation loss is estimated by Qsubl 5

The vapor transfer model developed by Thorpe and Mason (1966) and later modified by De´ry et al. (1998) is used to estimate daily sublimation losses as a function of snow properties and atmospheric conditions. Sublimation estimates in forested regions of the watershed are made only in the canopy, because we assume the effect of atmospheric turbulence on the ground surface

(13)

   Qr Ls 1 2prs  CNNu Ta Rv 3 Ta dm   ,  5 Ls Ls Ta dt  1 1 Rv CNNu Ta Rv 3 Ta N sh Dei (14) where 2pr (m) is the area function of a snow particle, s is the water vapor deficit, where

(12)

3) SUBLIMATION

dm N , dt (z)

where Qsubl (kg m22 s21) is the sublimation rate for a column of blowing snow over a horizontal land surface, dm/dt is the change in mass of a blowing snow particle as a result of sublimation per second, and N(z) is the number of snow particles per unit volume (m23). The number of snow particles depends on the particle shape a and radius r. A mean a 5 5 and r 5 100 mm (Pomeroy and Male 1988) are used. When using an a 5 5 and a 10-m wind speed of 15 m s21, De´ry et al. (1998) suggest that N(z) 5 9.09 3 107 m23. Because of the lack of wind data, we maintain the assumption that the average winter wind speed at 10-m height is 15 m s21. Thorpe and Mason (1966) estimate the change in mass of a blowing snow particle by

1 0.0285RainLAI 3 LAI  0.007LAI2 ,

where Rain is the precipitation that falls as rain. Von Hoyningen-Huene (1983) found that this function is unstable at precipitation values greater than 18 mm, resulting in anomalous interception values. Therefore, daily precipitation is set to a maximum of 18 mm for this calculation; all other precipitation is considered as throughfall. The intercepted rain store is determined for each day by adding the intercepted rain to the store.

VOLUME 10

s5

(e  ei) , ei

(15)

where e and ei are the vapor pressure, and its value at saturation over ice is determined by ei 5 8.55  0.20Ta 1 0.0457Ta2

(16)

  RH ei. e5 100

(17)

and

Here Qr is the radiation transferred to the particle, Qr 5 pr2(1 2 ap)Q (Schmidt 1991), with ap the shortwave * particle albedo of 0.5 (Schmidt et al. 1998); C is the thermal conductivity of air (2.4 3 1022 W m21 K21), Ls is the latent heat of sublimation (2.838 3 106 J kg21), Rv is the gas constant for water vapor (461.5 J kg21 K21), and D is the molecular diffusivity of water vapor in air

DECEMBER 2009

1439

MACDONALD ET AL.

(2.25 3 1025 m2 s21). The Nusselt number NNu and Sherwood number Nsh are N Nu and

0.5 N sh 5 1.79 1 (0.606N rs ),

(18)

where the Reynolds number Nre 5 (2rVr/V); Vr is the terminal ventilation velocity, where horizontal particle components are assumed equal to the horizontal wind speed (Schmidt 1982; De´ry and Taylor 1996), and V is the kinematic viscosity of air (1.53 3 1025 m2 s21; De´ry and Yau 1999).

4) SNOWMELT The snowmelt routine was adopted from Quick and Pipes (1977). For snowmelt to occur, the snowpack must ripen, which is determined by the ability of the snowpack to store cold. The snowpack cold storage (TREQ) is determined by TREQi 5 (MLTF 3 TREQi1 ) 1 T i ,

(19)

where MLTF is a decay constant (set to 0.85). When TREQ reaches 0 (enough energy is absorbed and the snowpack is ripe), melt can occur. Daily snowmelt (M; mm) values are calculated as a function of air temperature: M 5 PTM[Tmax 1 (TCEADJ 3 Tmin)],

(20)

where PTM is a point melt factor (mm day21 8C21) and TCEADJ is an energy partition multiplier [Eq. (21)], Tmax is maximum daily temperature, and Tmin is minimum daily temperature, TCEADJ 5

(Tmin 1 T) . (18 1 T)

TABLE 2. Comparison of date of complete melt at different PTM. Day (day of year) of complete melt Year

Observed

Simulated (PTM 1.0)

Simulated (PTM 1.8)

1977 1978 1979 1980 1981 1982 1983 1984 1985 1986

134 139 121 117 142 127 120 127 105 114

146 145 127 86 144 120 120 126 104 110

125 133 117 86 135 107 97 117 81 96

Q pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi * PET 5 0.038Q T 1 9.5  2.4 * Ra   RH , 1 0.075(T 1 20) 1  100

!2

(22)

where Ra is extraterrestrial radiation, determined as a function of Julian date. This equation is used because it does not require wind data and has been shown to provide accurate estimates of evapotranspiration when compared with the standardized Food and Agricultural Organization of the United Nations irrigation and drainage paper number 56 (FAO-56; Allen et al. 1998) Penman–Monteith scheme using a global climatic dataset (Valiantzas 2006). Actual evapotranspiration is restricted by the soil moisture; potential evapotranspiration is reduced when soil water declines below half of field capacity:

(21)

Wyman (1995) suggests a PTM of 1.8 for the Canadian Rocky Mountains. While maintaining the lapse rates of Pigeon and Jiskoot (2008), the timing and rate of melt at the Many Glacier SNOTEL site are used to calibrate PTM. On the basis of these results, a PTM of 1.0 appears to be more suitable to the study area. The rate of melt did not change significantly with changes in PTM. However, the timing of complete melt showed that a PTM of 1.0 provides the most suitable simulation of the date of complete melt, especially in years when complete melt occurred earlier (Table 2).

5) EVAPOTRANSPIRATION Daily potential evapotranspiration (PET) estimates are made only when the snowpack is depleted. When the snowpack remains, sublimation is calculated. A version of the Penman–Monteith evapotranspiration equation developed by Valiantzas (2006) is adopted:

 XK 5 2

1.5 SM , Soilmax

(23)

where XK is the water supply control on ET, SM is daily soil moisture (mm), and SoilMax is the soil field capacity (mm) for a particular TC. When SM is less than half of the storage capacity, ET is determined by ET 5 ETP 3 XK.

(24)

6) SOIL MOISTURE AND RUNOFF A daily soil water budgeting routine was developed for GENESYS. The soil water routine includes estimates of PET, ET, soil water storage, and controls on the ET/PET ratio where water supply limits ET. Estimates of the spatial variation in soil water processes were enabled using the Soil Survey Geographic (SSURGO) soils database (NRCS 2008).

1440

JOURNAL OF HYDROMETEOROLOGY

To account for changes in soil moisture, soils data for each TC are required. Soils data were only available for the eastern portion of the St. Mary River watershed. Therefore, land cover was used as a surrogate for soils in the western portion of the watershed. Using the same land cover grid as was used for the TC delineation, relationships between land cover type and soil type were used to extrapolate soil depth and water holding capacity values from the eastern portion of the watershed to the entire basin. Mean soil depth and field capacity values from each soil type were used for the analysis. A GIS overlay analysis—including land cover, mean soil depth, and mean field capacity grids—was used to extract mean soil depth and field capacity for each land cover type. Soils data were sparse in the eastern portion of the watershed above an elevation of 2700 m and a slope of 458. An analysis of land cover data showed that, although present, the density of vegetative cover is also low above this elevation and slope. It was, therefore, assumed that elevations above 2700 m and steeper than 458 slopes do not have significant water storage capacity. Therefore, the constraints ‘‘above 2700 m,’’ and ‘‘more than 458’’ were included in the analysis and assigned values of the barren rock–soil land cover type. For each land cover type, mean soil depth and field capacity values were determined. Using the land cover type as a surrogate for soil depth and field capacity, values from the eastern portion of the watershed were extrapolated to the western portion of the watershed. This is achieved by applying a soil depth and field capacity value to each land cover class, thereby representing the entire watershed. For each TC mean soil depth and field capacity values were determined. This enables the spatial representation of maximum soil field capacity over the entire drainage basin. Field capacity values over the watershed range from 0.0 to 199.4 mm, accounting for soil depth. The values associated with each TC are then used in the model to estimate changes to daily soil water loss through ET and gain from snowmelt and rain infiltration (IF). On days when either snowmelt or rainfall events occur and soil field capacity is exceeded, R is produced.

4. Model application The simulation was conducted from 1960 to 2001. This time frame was selected to incorporate the PRISM dataset and to enable application of future climate change scenarios. To evaluate the sensitivity of snowpack simulations to temperature and precipitation assumptions and assess the suitability of the GENESYS model in simulating daily snowpack, a comparison was made between

VOLUME 10

TABLE 3. Dates of selected MODIS imagery and percentage cloud cover. Date

Cloud cover (%)

16 Oct 2000 17 Nov 2000 17 Jan 2001 14 Mar 2001 15 Apr 2001 17 May 2001 12 Jul 2001 14 Sep 2001

0.3 0.3 5.8 1.1 0.2 1.9 2.5 0.0

the Many Glacier SNOTEL site and simulated values at the TC representing the SNOTEL site. These comparisons were made from the inception of the SNOTEL site in 1976 to the end of the 2001 snow year. To test how well the model simulates the spatial distribution of snow cover, a comparison was made with the MODIS/Aqua snow cover 8-day global 500-m grid, version 5 (Hall et al. 2007), for the dates presented in Table 3. The 8-day global 500-m grid provides one image for every eight days; however, not all data could be used. The dates shown in Table 3 were selected based on the condition that they cover both snow accumulation and ablation periods and that they contained less than 6% cloud cover. The spatial resolution of MODIS is 500 m 3 500 m. Therefore, to compare GENESYS snow cover to the MODIS data, snow cover surfaces were resampled to a spatial resolution of 500 m 3 500 m. Both MODIS and GENESYS snow surfaces were classified using a binary classification, where pixels with snow had a value of one and pixels without snow had a value of zero. For MODIS, pixels containing ice during the winter period were also included as snow-covered pixels. To evaluate the model performance, a similar method to that used by Garen and Marks (2005) was applied, where the percentage of pixels that were in agreement between MODIS and GENESYS snow surfaces for each of the eight dates was calculated. The Hanssen–Kuipers skill score (KSS) was also applied, which uses a contingency table where A (snow–snow), B (snow–no snow), C (no snow–snow), and D (no snow–no snow) classified pixels are evaluated. A perfect score for the test is a value of one. Because of the disproportionate number of pixels that are snow covered, we applied a form of the original KSS equation that enables equalization of snow and no-snow classified pixels (Woodcock 1976): KSS 5

A D 1 1. A1B C1D

(25)

DECEMBER 2009

MACDONALD ET AL.

TABLE 4. Average water balance for three elevation bands from 1961 to 2000.

Elevation

R (mm)

P (mm)

ET (mm)

SnowInt (mm)

RainInt (mm)

1500 2000 2500

690 1105 1560

1283 1661 1976

277 232 113

196 222 230

120 102 73

5. Results All components of the hydrological balance described in Eqs. (4) and (5) are shown for three elevation bands in Table 4. The simulation used for this demonstration of the hydrological balance applies the precipitation estimation method B, with temperature lapse rates derived by Pigeon and Jiskoot (2008). A 40-yr average for each of the hydrological balance components from 1961 to 2000 is presented, where elevation is in meters above sea level, runoff is the total runoff from rain and snowmelt, and P is both rain and snow. SnowInt is intercepted snow lost to sublimation. This is based on the assumption that intercepted snow that remains in the canopy is lost to sublimation, which has been observed by Pomeroy et al. (1998). RainInt is the amount of intercepted rain, ET is the water loss through evapotranspiration from the soil. Canopy rain loss from ET is accounted for in rainInt. Table 4 demonstrates that the structure of this model provides physically plausible estimates of hydrometeorological variables. It also shows the dependency of hydrometeorological estimates on elevation. It is shown that as elevation increases, there is an increase in the amount of precipitation received, a decrease in evapotranspiration, an increase in snow interception and sublimation, a decrease in rain interception, and an increase in runoff. To test precipitation estimate methods, simulated daily SWE was compared with observed daily SWE at

1441

the Many Glacier SNOTEL site at Many Glacier from 1 October 1976 to 26 April 2001 (Figs. 8–10). These dates correspond with the start and the end of the water years respectively. Figure 8 demonstrates that applying precipitation estimation method A results in relatively poor agreement between observed and simulated SWE values [r2 5 0.44, p , 0.0001, root-mean-square error (RMSE) 5 158 mm]. It is also evident that the magnitude of maximum SWE estimates are incorrect, with low SWE years being underestimated and some high SWE years being drastically overestimated (Fig. 8). Figure 9 is a comparison between observed and simulated SWE at Many Glacier using precipitation method B. There is good agreement between observed and simulated SWE values (r2 5 0.72, p , 0.0001, RMSE 5 88 mm). This result is a significant improvement relative to SWE estimates using precipitation method A. Figure 10 shows that using precipitation method C, where PRISM data are applied, results in the best agreement between observed and simulated SWE values at Many Glacier (r2 5 0.73, p , 0.0001, RMSE 5 73 mm). Simulated SWE values using PRISM input compare well with simulated SWE values using precipitation method B, showing that accounting for differences between the precipitation climatology of mountainous and transitional prairie zones is important. The sensitivity of temperature lapse rates on daily SWE simulations was tested for the 1982 snow year (Fig. 11). This year was selected because it was the most accurately simulated over the time series (r2 5 0.98, p , 0.0001). Linear regression was applied to test the difference between observed SWE and simulated SWE using all three lapse rates. The coefficient of determination did not differ between simulations. The RMSE did, however, differ slightly with both of the static lapse rates having a RMSE of 30 mm, whereas the daily lapse rates resulted in a RMSE of 27 mm.

FIG. 8. Observed vs simulated SWE using precipitation method A.

1442

JOURNAL OF HYDROMETEOROLOGY

VOLUME 10

FIG. 9. Same as Fig. 8, but using precipitation method B.

The spatial snow cover extent simulated by GENESYS compared well with MODIS snow cover for the dates selected. The percentage of correctly classified pixels ranged from 73% (17 January 2001) to 100% (12 July 2001 and 14 September 2001), with an average of 90% correctly classified pixels for all eight dates. The KSSs ranged from 0.16 (17 January 2001) to 1.0 (12 July 2001 and 14 September 2001; Table 5). Figure 12 compares MODIS snow cover with simulated snow cover for three dates. The dates were selected based on their representativeness of the snow accumulation and ablation periods.

6. Discussion This study has demonstrated the applicability of the GENESYS model in estimating snowpack over a mountainous watershed. This physically based model provides an alternative to spatial interpolation and can operate

at a fine spatial scale with a relatively low level of required input data. The spatial structure of the model enables applicability to a variety of landscapes and provides a useful tool for spatial hydrometeorological simulations in watersheds with little observed data. The hydrological balance simulated by GENESYS is physically realistic. With an average of 15% over three elevation bands, the proportion of the total annual precipitation lost to sublimation in the canopy compares well with Strasser et al. (2007), Hood et al. (1999), and Zhang et al. (2004). The intercepted rain loss accounts for an average of 6% of the annual precipitation received over the three elevation bands. This seems reasonable given the relatively large proportion of the annual precipitation received as snow. The decrease in evapotranspiration with elevation, despite an increase in precipitation, demonstrates the sensitivity of the evapotranspiration routine to temperature estimates. This is expected given temperature is a key variable in determining

FIG. 10. Same as Fig. 8, but using precipitation method C.

DECEMBER 2009

1443

MACDONALD ET AL.

FIG. 11. Comparison of the effect of lapse rate on SWE simulations.

vapor pressure deficit (Valiantzas 2006). Table 4 shows that runoff estimates are largely a function of the precipitation received. Table 4 also demonstrates that the runoff contributions to streamflow from higher elevations in the St. Mary River watershed are likely important because a greater proportion of the annual precipitation received is accounted for by runoff. With realistic simulations of the average annual hydrological balance, we have confidence that the GENESYS model can account for snow accumulation and ablation. A detailed analysis of snowpack simulations provides further insight into the model’s ability to simulate these processes at a finer time step. The method used to predict precipitation over the watershed can provide realistic monthly precipitation estimates, which can be used to simulate SWE at a daily time step. It is plausible that this method can be applied to other watersheds with low-elevation climate stations and relatively few snowpack observations for large-scale hydrological simulations. The PRISM 1971–2000 monthlymean precipitation dataset was shown to be a reliable source for spatial precipitation inputs to snowpack simulations. The key advantage of these two methods of precipitation estimation is that they are able to describe the climatic characteristics of mountainous regions of the watershed. This study, and others (Bales et al. 2006), has shown that mountainous regions experience very different precipitation regimes relative to low-elevation mountain– prairie transitional zones. Figure 8 demonstrates that by not accounting for this seasonal difference, the method used by LBJLK fails to represent precipitation at the Many Glacier site. Figures 9 and 10 show that accounting for seasonal differences through the method developed in this paper and the PRISM method enable more representative simulations of daily SWE. The interannual variability in the accuracy of SWE simulations demonstrates that the meteorological conditions at the St. Mary climate station are not always representative of the Many

Glacier SNOTEL site. This variability also suggests that applying a mean precipitation–elevation relationship does not account for year-to-year changes in the relationship. However, given the available data, the method provides reasonable estimates of SWE overall. Assumptions in estimating temperatures have little effect on the goodness of fit of the SWE simulations at Many Glacier (Fig. 11). However, the timing of snow accumulation and ablation are highly sensitive to the lapse rates used. Model runs that apply the lapse rates of LBJLK and Pigeon and Jiskoot (2008) overestimate the snow accumulation period and total SWE for the 1981 snow year. Figure 11 shows that there is little difference between these lapse rates in their effects on simulated daily SWE during the snow accumulation period. During the snowmelt period, however, temperature simulations that use the lapse rates derived by Pigeon and Jiskoot (2008) result in a more accurate simulation of the date of complete melt. The lapse rates derived using NCEP performed slightly better than the static lapse rates. The NCEP lapse rates also resulted in the most accurate simulation of the snow accumulation period and total SWE. The NCEP lapse rates did not, however, result in the most accurate simulation of the date of complete melt.

TABLE 5. Percent correctly classified snow pixel and KSS test results for comparison of simulated snow cover and MODIS snow cover maps for eight dates. Date

KSS

Correct (%)

16 Oct 2000 17 Nov 2000 17 Jan 2001 14 Mar 2001 15 Apr 2001 17 May 2001 12 Jul 2001 14 Sep 2001

0.36 0.96 0.15 0.99 0.22 0.60 1.0 1.0

78 95 73 98 97 80 100 100

1444

JOURNAL OF HYDROMETEOROLOGY

VOLUME 10

FIG. 12. Spatial comparison of MODIS snow cover images and simulated snow cover using the GENESYS models for three dates during the 2000 water year.

The snow accumulation and ablation periods were also assessed spatially by comparing MODIS snow cover data with simulated snow cover for eight dates. Snow cover simulations for 17 November 2000, 17 January 2001, and 17 May 2001 compare well with MODIS snow cover (Fig. 12). KSS values have more variability than the percentage correctly classified pixel calculations. This is likely because percentage of snow pixels correctly classified does not account for those pixels that are nonsnow covered and classified as snow or vice versa. This result also demonstrates the sensitivity of the KSS test. The 17 November 2000 simulation underestimates snow accumulation in the upper watershed, with only the highest elevations being covered by snow. This underestimate is likely a function of warm temperatures at the St. Mary climate station, which would result in less snow accumulation. This is supported by the 17 January simulation. Without snow cover, surface heating over St. Mary from a reduced albedo relative to a snowcovered surface likely results in an oversimulation of temperature at low elevations. The 17 May simulation has good visual agreement with MODIS, although only 80% of the pixels were correctly classified and the KSS value is 0.60. The results of this spatial comparison demonstrate the ability of the model to use very little input data to reasonably represent spatial snow cover over a mountainous watershed.

7. Limitations Independent of the precipitation estimation method, this study has demonstrated the difficulties in using lowelevation climate records for driving spatially distributed hydrological models in complex terrain. Using a single low-elevation station to determine a watershed scale hydrological balance is subject to significant influence from local hydrometeorological conditions. To mitigate this problem in the future, further meteorological monitoring is necessary. The dependence of precipitation on topography is also only partially accounted for when using elevation as the sole predictor. The precipitation–elevation relationship used in this study explains a maximum of 44% of the variability in mean winter precipitation. The model would likely explain more of the variability in precipitation if other topographical predictors could be used. For instance, Anderton et al. (2004), Basist et al. (1994), Daly et al. (1994, 2002, 2007, 2008), and Marquinez et al. (2003) have shown that using a combination of variables including slope, aspect, coastal proximity, and orientation can greatly increase the predictive power of regressionbased estimates of precipitation. Relationships between precipitation and slope–aspect were not significant for the Preston snow course (LBJLK), likely because the data do not represent enough variability in slope and aspect.

DECEMBER 2009

MACDONALD ET AL.

8. Summary Limited data can be used to simulate the spatial distribution of hydrometeorological characteristics of a mountainous watershed. The GENESYS model can provide physically realistic spatial estimates of snow cover and snowpack accumulation and ablation in a mountainous watershed. This study, however, also supports suggestions by Daly et al. (2007) and Bales et al. (2006) that the key limitation to modeling in mountainous regions is the lack of observed data for model calibration and validation. For simulating large watersheds, it is necessary to identify the critical components of large-scale snow hydrology and represent these processes with appropriate equations that are able to resolve the complexities at the relevant spatial and temporal scales (Woo and Marsh 2005). With increased monitoring, physical studies would be enabled, describing important processes such as wind redistribution of snow (e.g., Anderton et al. 2004) and cold air drainage (e.g., Lundquist and Cayan 2007). Further studies using the GENESYS model will attempt to quantify these micrometeorological processes and their subsequent effects on modeling at the watershed scale, thereby determining the level at which to monitor ecosystem processes for management applications in relatively large watersheds. Acknowledgments. Funding support from the Alberta Ingenuity Centre for Water Research and the Natural Sciences and Engineering Research Council of Canada is much appreciated. Dr. Dan Fagre, USGS, kindly provided the snow survey data used in this study. The advice of the reviewers is very much appreciated, as they have greatly improved the manuscript. REFERENCES Allen, R. G., L. S. Pereira, D. Raes, and M. Smith, 1998: Crop evapotranspiration: Guidelines for computing crop water requirements. Food and Agriculture Organization of the United Nations Irrigation and Drainage Paper 56, 328 pp. Anderton, S. P., S. M. White, and B. Alvera, 2004: Evaluation of spatial variability in snow water equivalent for a high mountain catchment. Hydrol. Processes, 18, 435–453. Bales, R. C., N. P. Molotch, T. H. Painter, M. D. Dettinger, R. Rice, and J. Dozier, 2006: Mountain hydrology of the western United States. Water Resour. Res., 42, 13. Band, L. E., D. L. Peterson, S. W. Running, J. Coughlan, R. Lammers, J. Dungan, and R. Nemani, 1991: Forest ecosystem processes at the watershed scale: Basis for distributed simulation. Ecol. Modell., 56, 171–196. ——, P. Patterson, R. Nemani, and S. W. Running, 1993: Forest ecosystem processes at the watershed scale: Incorporating hillslope hydrology. Agric. For. Meteor., 63, 93–126. Barnett, T. P., J. C. Adam, and D. P. Lettenmaier, 2005: Potential impacts of a warming climate on water availability in snowdominated regions. Nature, 438, 1–7.

1445

Basist, A., G. D. Bell, and V. Meentemeyer, 1994: Statistical relationships between topography and precipitation patterns. J. Climate, 7, 1305–1315. Beniston, M., 2003: Climatic change in mountain regions: A review of possible impacts. Climatic Change, 59, 5–31. Daly, C., 2006: Guidelines for assessing the suitability of spatial climate data sets. Int. J. Climatol., 26, 707–721. ——, R. P. Neilson, and D. L. Phillips, 1994: A statistical–topographic model for mapping climatological precipitation over mountainous terrain. J. Appl. Meteor., 33, 140–158. ——, W. P. Gibson, G. H. Taylor, G. L. Johnson, and P. P. Pasteris, 2002: A knowledge-based approach to the statistical mapping of climate. Climate Res., 22, 99–113. ——, J. W. Smith, J. I. Smith, and R. B. McKane, 2007: Highresolution spatial modeling of daily weather elements for a catchment in the Oregon Cascade Mountains, United States. J. Appl. Meteor. Climatol., 46, 1565–1586. ——, M. Halbleib, J. I. Smith, W. P. Gibson, M. K. Doggett, G. H. Taylor, J. Curtis, and P. A. Pasteris, 2008: Physiographically sensitive mapping of climatological temperature and precipitation across the conterminous United States. Int. J. Climatol., 28, 2031–2064. De´ry, S. J., and P. A. Taylor, 1996: Some aspects of the interaction of blowing snow with the atmospheric boundary layer. Hydrol. Processes, 10, 1345–1358. ——, ——, and J. Xiao, 1998: The thermodynamic effects of sublimating, blowing snow in the atmospheric boundary layer. Bound.-Layer Meteor., 89, 251–283. ——, and M. K. Yau, 1999: A bulk blowing snow model. Bound.Layer Meteor., 93, 237–251. Diaz, H. F., 2005: Monitoring climate variability and change in the western United States. Global Change and Mountain Regions, U. M. Huber, H. K. M. Bugmann, and M. A. Reasoner, Eds., Advances in Global Change Research, Vol. 23, Springer, 636–656. Field, C. B., and Coauthors, 2007: North America. Climate Change 2007: Impacts, Adaptation and Vulnerability, M. L. Parry et al., Eds., Cambridge University Press, 617–652. Fontaine, T. A., T. S. Cruickshank, J. G. Arnold, and R. H. Hotchkiss, 2002: Development of a snowfall-snowmelt routine for mountainous terrain for the soil water assessment tool (SWAT). J. Hydrol., 262, 209–223. Garen, D. C., and D. Marks, 2005: Spatially distributed energy balance snowmelt modelling in a mountainous river basin: Estimation of meteorological inputs and verification of model results. J. Hydrol., 315, 126–153. Glassy, J. M., and S. W. Running, 1994: Validating diurnal climatology logic of the MT-CLIM model across a climatic gradient in Oregon. Ecol. Appl., 4, 248–257. Groisman, P., T. R. Karl, and T. W. Knight, 1994: Observed impact of snow cover on the heat balance and the rise of continental spring temperatures. Science, 263, 198–200. Hall, D. K., G. A. Riggs, and V. V. Salomonson, 2007: MODIS/ Aqua snow cover 8-day L3 global 500m grid, version 5. National Snow and Ice Data Center, Boulder, CO, digital media. [Available online at http://nsidc.org/data/myd10a2v5.html.] Hamlet, A. F., and D. P. Lettenmaier, 1999: Effects of climate change on hydrology and water resources in the Columbia basin. J. Amer. Water Resour. Assoc., 35, 1597–1623. Hanson, C. L., 1982: Distribution and stochastic generation of annual and monthly precipitation on a mountainous watershed in southwest Idaho. J. Water Resour. Assoc., 18, 875–883.

1446

JOURNAL OF HYDROMETEOROLOGY

Hedstrom, N. R., and J. W. Pomeroy, 1998: Measurements and modelling of snow interception in the boreal forest. Hydrol. Processes, 12, 1611–1625. Hood, E., M. Williams, and D. Cline, 1999: Sublimation from a seasonal snowpack at a continental, mid-latitude alpine site. Hydrol. Processes, 13, 1781–1797. Hungerford, R. D., R. R. Nemani, S. W. Running, and J. C. Coughlan, 1989: MTCLIM: A mountain microclimate simulation model. USDA Forest Service Research Paper INT-414, 52 pp. Kane, D. L., L. D. Hinzman, C. S. Benson, and G. E. Liston, 1991: Snow hydrology of a headwater Arctic basin 1. Physical measurements and process studies. Water Resour. Res., 27, 1099–1109. Kienzle, S. W., 2008: A new temperature based method to separate rain and snow. Hydrol. Processes, 22, 5067–5085. Lapp, S., J. Byrne, S. Kienzle, and I. Townshend, 2002: Linking global circulation model synoptics and precipitation for western North America. Int. J. Climatol., 22, 1807–1817. ——, J. M. Byrne, I. Townshend, and S. W. Kienzle, 2005: Climate warming impacts on snowpack accumulation in an alpine watershed. Int. J. Climatol., 25, 521–526. Leavesley, G. H., R. W. Lichty, B. M. Troutman, and L. G. Saindon, 1983: Precipitation- runoff modeling system: User’s manual. USGS Water-Resources Investigations Rep. 83-4238, 207 pp. Lehning, M., I. Volksch, D. Gustafsson, T. A. Nguyen, M. Stahli, and M. Zappa, 2006: ALPINE3D: A detailed model of mountain surface processes and its application to snow hydrology. Hydrol. Processes, 20, 2111–2128. Leung, L. R., and M. S. Wigmosta, 1999: Potential climate change impacts on mountain watersheds in the Pacific northwest. J. Amer. Water Resour. Assoc., 35, 1463–1471. Liston, G. E., and K. Elder, 2006a: A distributed snow-evolution modeling system (SnowModel). J. Hydrometeor., 7, 1259–1276. ——, and ——, 2006b: A meteorological distribution system for high-resolution terrestrial modeling (MicroMet). J. Hydrometeor., 7, 217–234. Lundquist, J. D., and D. R. Cayan, 2007: Surface temperature patterns in complex terrain: Daily variations and long-term change in the central Sierra Nevada, California. J. Geophys. Res., 124, D11124, doi:10.1029/2006JD007561. Marks, D., J. Dozier, and R. E. Davis, 1992: Climate and energy exchange at the snow surface in the alpine region of the Sierra Nevada 1. Meteorological measurements and monitoring. Water Resour. Res., 28, 3029–3042. Marquinez, J., J. Lastra, and P. Garcia, 2003: Estimation models for precipitation in mountainous regions: The use of GIS and multivariate analysis. J. Hydrol., 270, 1–11. Mather, J. R., and G. A. Yoshioka, 1968: The role of climate in the distribution of vegetation. Ann. Assoc. Amer. Geogr., 58, 29–41. McKenzie, D., D. W. Peterson, D. L. Peterson, and P. E. Thornton, 2003: Climatic and biophysical controls on conifer species distributions in mountain forests of Washington State, USA. J. Biogeogr., 30, 1093–1108. Mote, P. W., A. F. Hamlet, M. P. Clark, and D. P. Lettenmaier, 2005: Declining mountain snowpack in western North America. Bull. Amer. Meteor. Soc., 86, 39–49. NCDC, cited 2006: Web climate service. [Available online at www.ncdc.noaa.gov/oa/climate/stationlocator.html.] NCEP, cited 2008: NCEP Global Ocean Data Assimilation System (GODAS). [Available online at http://www.cdc.noaa.gov/.]

VOLUME 10

NRCS, cited 2006: SNOTEL data and products. [Available online at http://www.wcc.nrcs.usda.gov/snow/.] ——, cited 2008: SSURGO data and products. [Available online at http://soildatamart.nrcs.usda.gov/.] Pigeon, K., and H. Jiskoot, 2008: Meteorological controls on snowpack formation and dynamics in the southern Canadian Rocky Mountains. Arct. Antarct. Alp. Res., 40, 716–730. Pomeroy, J. W., and D. H. Male, 1988: Optical properties of snow. J. Glaciol., 34, 3–9. ——, J. P. Parviainen, N. Hedstrom, and D. M. Gray, 1998: Coupled modelling of forest snow interception and sublimation. Hydrol. Processes, 12, 2317–2337. Quick, M. C., and A. Pipes, 1977: UBC Watershed Model. Hydrol. Sci. Bull., 306, 215–233. Roy, D. P., J. S. Devadiga, R. E. Wolfe, M. Zheng, and J. Descloitres, 2002: The MODIS land product quality assessment approach. Remote Sens. Environ., 83, 62–76. Running, S. C., R. R. Nemani, and R. D. Hungerford, 1987: Extrapolation of synoptic meteorological data in mountainous terrain and its use for simulating forest evapotranspiration and photosynthesis. Can. J. For. Res., 17, 472–483. Schmidt, R. A., 1982: Vertical profiles of wind speed, snow concentrations, and humidity in blowing snow. Bound.-Layer Meteor., 23, 223–246. ——, 1991: The sublimation of snow intercepted by an artificial conifer. Agric. For. Meteor., 54, 1–27. ——, C. A. Troendle, and J. R. Meiman, 1998: Sublimation of snowpacks in subalpine conifer forests. Can. J. For. Res., 28, 501–513. Selkowitz, D. J., D. B. Fagre, and B. A. Reardon, 2002: Interannual variations in snowpack in the Crown of the Continent Ecosystem. Hydrol. Processes, 54, 3651–3665. Stephenson, N. L., 1990: Climatic control of vegetation distribution: The role of the water balance. Amer. Nat., 135, 649–670. Strasser, U., M. Bernhardt, M. Weber, G. E. Liston, and W. Mauser, 2007: Is snow sublimation important in the alpine water balance? Cryosphere Discuss., 1, 303–350. Thorpe, A. D., and B. J. Mason, 1966: The evaporation of ice spheres and ice crystals. Br. J. Appl. Phys., 17, 541–548. USGS, cited 2006: Metadata for national land cover data for Montana. [Available online at http://nris.mt.gov/nsdi/nris/nlcdgrid.html.] Valiantzas, J. D., 2006: Simplified versions for the Penman evaporation equation using routine weather data. J. Hydrol., 331, 690–702. Von Hoyningen-Huene, J., 1983: Die Interzeption des Niederschlages in landwirtschaftlichen Pflanzenbesta¨nden. Dtsch. Verb. Wasserwirtsch. Kulturbau, 57, 1–66. Woo, M. K., and P. Marsh, 2005: Snow, frozen soils and permafrost hydrology in Canada, 1999-2002. Hydrol. Processes, 19, 215–229. Woodcock, F., 1976: The evaluation of yes/no forecasts for scientific and administrative purposes. Mon. Wea. Rev., 104, 1209–1214. Wyman, R. R., 1995: Modeling snowpack accumulation and depletion. Mountain Hydrology: Peaks and Valleys in Research and Applications, B. T. Guy and J. Barnard, Eds., Canadian Water Resources Association, 23–30. Zhang, T., R. L. Armstrong, and J. Smith, 2003: Investigation of the near-surface soil freeze-thaw cycle in the contiguous United States: Algorithm development and validation. J. Geophys. Res., 108, 8860, doi:10.1029/2003JD003530. Zhang, Y., K. Suzuki, T. Kadota, and T. Ohata, 2004: Sublimation from snow surface in southern mountain taiga of eastern Siberia. J. Geophys. Res., 109, D21103, doi:10.1029/2003JD003779.