Using ground penetrating radar to assess the variability of snow water equivalent and melt in a mixed canopy forest, Northern Colorado

Abstract

Snow is an important environmental variable in headwater systems that controls hydrological processes such as streamflow, groundwater recharge, and evapotranspiration. These processes will be affected by both the amount of snow available for melt and the rate at which it melts. Snow water equivalent (SWE) and snowmelt are known to vary within complex subalpine terrain due to terrain and canopy influences. This study assesses this variability during the melt season using ground penetrating radar to survey multiple plots in northwestern Colorado near a snow telemetry (SNOTEL) station. The plots include south aspect and flat aspect slopes with open, coniferous (subalpine fir,Abies lasiocarpa and engelman spruce, Picea engelmanii), and deciduous (aspen, populous tremuooides) canopy cover. Results show the high variability for both SWE and loss of SWE during spring snowmelt in 2014. The coefficient of variation for SWE tended to increase with time during snowmelt whereas loss of SWE remained similar. Correlation lengths for SWE were between two and five meters with melt having correlation lengths between two and four meters. The SNOTEL station regularly measured higher SWE values relative to the survey plots but was able to reasonably capture the overall mean loss of SWE during melt. Ground Penetrating Radar methods can improve future investigations with the advantage of non-destructive sampling and the ability to estimate depth, density, and SWE.

Keywords

headwaters / snowmelt / snow water equivalent / ground penetrating rdar / SNOTEL

Introduction

Snow is an important part of the hydrologic cycle and ecosystem dynamics for headwater systems. Many headwater regions in mountainous terrain have snowmelt dominated hydrographs that provide valuable in-stream water resources (Daly et al., 2000; Bales et al., 2006; Seyfried et al., 2009; Rice et al., 2011), pose a risk of flood damage (Graybeal and Leathers, 2006; Zhao et al., 2009; Fang et al., 2015), and recharge groundwater resources (Flint et al., 2008; Clilverd et al., 2011; Cao et al., 2012). Snow processes such as accumulation and melt will vary with elevation, topographic influences, and canopy effects (Elder et al., 1991; Blöschl and Kirnbauer, 1992; Molotch et al., 2005; Fassnacht and Derry, 2010; Clark et al., 2011; López-Moreno et al., 2011, 2013; Richer et al., 2013; Molotch and Meromy, 2014; Sexstone and Fassnacht, 2014; Fassnacht and Hultstrand, 2015). These factors and the resulting snowmelt variability have proven important for regular diurnal and seasonal fluctuations in streamflow (Lundquist and Cayan, 2002; Lundquist et al., 2005; Jencso and McGlynn, 2011; Mutzner et al., 2015) and stream connectivity to the surrounding landscape (Liu et al., 2004; McNamara et al., 2005; Seyfried et al., 2009; Eiriksson et al., 2013). Furthermore, headwater systems have shown variable trends in snow accumulation and melt patterns as a result of climate change (Adam et al., 2009; Clow, 2010; Harpold et al., 2012; Fassnacht and Hultstrand, 2015; Fassnacht et al., 2016). With climate change producing variable patterns and trends it is important to understand the underlying processes that drive hydrological dynamics in a headwater system during spring snowmelt.

Spring snow surveys have shown the dynamic spatio-temporal variability of snow depth and other properties such as density and the resulting snow water equivalent (SWE) in complex terrain (Elder et al., 1991; Cline et al., 2009; Elder et al., 2009; López-Moreno et al., 2011, López-Moreno et al., 2013; Sexstone and Fassnacht, 2014). Terrain-driven variability in these properties occurs from the centimeter to the kilometer scale (Blöschl, 1999; Fassnacht et al., 2009; Richer et al., 2013; Sexstone and Fassnacht, 2014) with important implications towards monitoring snow and snowmelt using single-point measurements for streamflow prediction purposes (López-Moreno et al., 2011; López-Moreno et al., 2013) and for soil moisture or groundwater recharge estimates (Harpold et al., 2015; Webb et al., 2015). The variability of snowmelt has shown standard point measurements such as snow pillows can be inaccurate, particularly early in the spring snowmelt season (Johnson and Schaefer, 2002). Additionally, the spatial variability of depth during spring snowmelt can influence the preferential flow of snowmelt water, altering the albedo of the snowpack at correlation lengths shown to be between five and seven meters in alpine environments (Sommerfeld et al., 1994; Williams et al., 1999) but correlation lengths have not been investigated in subalpine environments with additional complexity from tree canopies.

Tree canopy influences have been shown to cause higher variability of snow depth and snowmelt (Storck et al., 2002; Winkler et al., 2005; López-Moreno and Latron, 2008; Andreadis et al., 2009; Varhola et al., 2010; Broxton et al., 2015; Moeser et al., 2016) with implications towards infiltration patterns (Webb et al., 2015; Lundberg et al., 2016). This is largely due to the shortwave radiation shading that occurs as a result of coniferous canopies that remain through the winter season (Musselman et al., 2008, 2012; Molotch et al., 2009). Trees without canopies in winter due to a deciduous nature or a result of forest fire have also displayed effects on snow through wind resistance in addition to shortwave and longwave radiation changes in the energy balance (Ebel et al., 2012; Harpold et al., 2014). These effects result in spatial variability of snow ablation in subalpine terrain that will affect stream connectivity to landscapes (McNamara et al., 2005; Jencso et al., 2009; Seyfried et al., 2009), spatial variability of groundwater recharge (Clilverd et al., 2011; Magnusson et al., 2014; Webb et al., 2015), and soil moisture dynamics important to plant production (Williams et al., 2009a; Williams et al., 2009b; Bales et al., 2011; Harpold et al., 2015). Quantification of these hydrological processes will be affected by estimates of snowmelt variability.

Many methods to assess the spatio-temporal variability of snowmelt have been utilized over the years. Recently, Ground Penetrating Radar (GPR) methods have been developed as a non-destructive method to measure snow depth (Gusmeroli and Grosse, 2012), density (Previati et al., 2011), wetness (Granlund et al., 2009; Mitterer et al., 2011), and even interpreted for stratigraphy information (Heilig et al., 2009). GPR data offers the benefit of non-destructive methods to obtain near continuous spatial measurements of depth and estimates of SWE in a relatively short period of time. The aim of this paper is to utilize GPR methods to evaluate patterns of SWE and changes in SWE during spring snowmelt in subalpine terrain. The dataset collected to accomplish this utilizes GPR technology for spatial surveys along transects that cross over multiple land cover types near a Snow Telemetry (SNOTEL) site that offers continuous snow and meteorological observations. These surveys were repeated during spring snowmelt and the resulting data analyzed to: 1) determine how representative single point and/or plot measurements of SWE and changes in SWE are for similar plot types nearby during melt, 2) compare the loss of SWE during snowmelt recorded at the SNOTEL site to the surrounding areas of variable slope aspect and vegetative cover, and 3) determine the range of correlation in SWE and changes in SWE for varying land cover types.

Methods

Site description

The location of this investigation is the Dry Lake (DL) study site that is approximately 6.5 km northeast of Steamboat Springs, Colorado located in Routt National Forest; data collection was within an area of interest less than 0.1 km² in area (Fig. 1). The elevation of GPR surveys range from 2500 m to 2540 m with slope angles from 1° to 15° as determined from a 10 m digital elevation model (USGS, 2015). The site has a mix of deciduous (aspen, Populous tremuloides) and evergreen forest (subalpine fir, Abies lasiocarpa; engelman spruce, Picea engelmanii) with a majority of the vegetation growing near a small stream. Slope aspects and canopy covers surveyed were classified as: south aspect with open canopy (SO), flat aspect with open canopy (FO), flat aspect with deciduous canopy (FD), and flat aspect with coniferous canopy (FC).

At approximately 2540 m elevation along an exposed ridge at the top of the south aspect slope is a Remote Automated Weather Station (RAWS). For this study, the RAWS location is only used for fixed landmarks to assist with survey transects as the station does not record any snow data. However, the Dry Lake SNOTEL station is located approximately 120 m to the south-southwest of the RAWS at a lower elevation of 2510 m in an open area near deciduous trees (Fig. 1). Data from the SNOTEL station include daily observations of air temperature, precipitation, relative humidity, snow depth, SWE, soil moisture and temperature, solar radiation, and wind speed and direction. All data are available from 2003 to present, with SWE and precipitation being collected since 1979. The SNOTEL data show peak SWE occurs April 5 on average with a 35 year median peak of 570 mm and a mean of 590 mm.

Data collection and analysis

Four Transects were established using stationary positions such as the southeast corner of the SNOTEL station, trees marked with flagging, and the RAWS site. The use of a differential global positioning system (GPS) was not available for this study, thus a handheld GPS unit was used to approximate locations of transects for mapping purposes. The use of fixed landmarks and flagging offer locations with accuracy at the tens of centimeter scale relative to the handheld GPS accuracy at the meters scale. Landmarks were established within line of site from previous landmarks along transects to ensure proper direction of travel. Maintaining direction of travel at this location had few difficulties with skis over the snow and minimal distances between landmarks. Total transect lengths ranged from almost 120 m up to 200 m. Each transect was established on April 20, 2014 when first surveyed with the following surveys occurring in two week intervals (May 4 and May 18). The GPR data was analyzed using the two way travel time (TWT) of the radar wave from the reflected ground surface (Lundberg et al., 2006; Heilig et al., 2009; Mitterer et al., 2011). The wave speed will depend upon the dielectric properties of the medium using:

(1)v={\displaystyle \frac{c}{\sqrt{{\epsilon }_r}}},

where v is the velocity of the electromagnetic wave (m/ns), c is the speed of light in a vacuum (0.3 m/ns), and e_r is the relative dielectric permittivity of the medium. For snow, e_r will depend upon the dry density of the snowpack and the liquid water content because liquid water and ice have different permittivities (e.g.,Heilig et al., 2015). Depth can then be estimated from an estimate of the velocity using:

(2)d={\displaystyle \frac{TWT\cdot v}{\mathrm{2}}},

where d is the snow depth in meters. For this study, known snow depths at three locations were used to calibrate the GPR wave speed estimates. Wave speeds in this study were adjusted within a range from a minimum allowable velocity of 0.13 m/ns to a maximum of 0.24 m/ns (Heilig et al., 2009; Mitterer et al., 2011). One depth measurement utilizes the SNOTEL station data for FO conditions whereas the other two were from snow pits dug along transect 1 under FD conditions and transect 2 under SO conditions (Fig. 1). SWE was then estimated for transects from bulk density measurements at these snow pits and SNOTEL station by multiplying depth by the appropriate density associated with the land cover and slope aspect type. SWE estimates across each plot from GPR interpretation assumes that density will not vary significantly for each slope and land cover type (López-Moreno et al., 2013) and results were compared to seven additional density measurements collected during the same survey days within 500 m of the transects as part of a separate study (Webb and Fassnacht, 2016). GPR surveys were conducted early in the morning, prior to the onset of the days melt to reduce errors from preferential flow of liquid water content in the snowpack.

The GPR unit used was the pulse EKKO pro with 1.0 GHz antennas, manufactured by Sensors & Software, placed in a sled for travel across the snow surface. This frequency of radar is the highest recommended to penetrate snowpacks up to two meters deep that may contain liquid water (Heilig et al., 2009; Mitterer et al., 2011; Singh et al., 2011). The weight of the sled caused it to compress the surface of the snow by five centimeters. Because of this compression, all depth estimates were adjusted to account for this. Fiducial points were recorded at the established landmarks on each transect for aligning datasets during analysis. Step sizes for measurements were 10 cm with six times stacking of recordings. Distances along all transects were measured using a “big wheel odometer” from Sensors & Software that was calibrated along a 50 m distance of snow surface upon approach to the study site. Visualization of survey data was accomplished using EKKOView Deluxe software. Filtering of the data for interpretation was done using the dewow function that removes the long waved part of the signal caused by electromagnetic induction followed by spherical and elliptical correction (SEC) to compensate for losses and dissipation of energy. SWE estimates were confirmed from SNOTEL observations, snow pit measurements at two locations directly adjacent to transects, seven more snow pits with variable canopy cover within 500 m as part of a separate study (Webb and Fassnacht, 2016), and radar wave velocity estimates from hyperbolic diffraction curvature when diffractions were present at the ground surface.

Interpreted data was used to create SWE profiles along each transect at a measurement spacing of 10 cm and interpreted for SWE with 1.0 mm precision. Each transect was then segregated into datasets for each individual survey plot. Loss of SWE profiles were calculated for melt period 1 (May 4 SWE profile subtracted from April 20 profile) and melt period 2 (May 18 SWE subtracted from May 4). This creates two datasets for statistical analysis, SWE and loss of SWE. Histograms of the datasets were created for each land cover type in addition to variograms to estimate correlation lengths of SWE and loss of SWE for each individual survey plot. The loss of SWE for varying land cover type and slope aspect was then further compared to the measured losses at the SNOTEL station to assess how well it captures melt characteristics of the surrounding snowpack.

Results

GPR surveys resulted in clear reflections from the ground surface beneath the snowpack at the DL study site for all three survey periods (Fig. 2, example shown: May 4, transect 1). For each survey date, known snow depths and densities at three locations were used to calibrate GPR wave velocities for FO, FD, and SO conditions; the FC survey plot was assumed to have the same wave velocity and snow density as the adjacent FD plot. The minimum calibrated wave velocity (0.13 m/ns) occurred for SO conditions and the maximum wave velocity (0.24 m/ns) occurred under FD conditions. These represent the minimum and maximum velocities allowed in the calibration method. However, the majority of the calibrated wave velocities were between 0.15 m/ns and 0.23 m/ns for the entire observation period. The low wave velocity on the SO slope indicates higher liquid water content and/or higher relative density for the snowpack whereas the higher wave velocities indicate a snowpack that is less dense and/or has a lower liquid water content. The higher velocities occur in plot types with lower bulk density observations from snow pits and low velocities in higher bulk density snow pit observations (Fig. 3). The density observations used for calibration compared well with other snow pits in the area and differences were generally less than 10% and maximum differences for SO pits at 10%, FO compared to SNOTEL at 15%, and FD compared to FC at 13% (Fig. 3).

The resulting SWE profiles display a maximum observed SWE of 825 mm, 634 mm, and 378 mm for April 20, May 4, and May 18, respectively; all occurring in FO survey plots (Fig. 4). The minimum measured SWE was 159 mm, 0 mm and 0 mm for April 20, May 4 and May 18, respectively; occurring in SO survey plots (Fig. 4). The loss of SWE for the surveys resulted in maximum losses of 444 mm and 480 mm for melt periods 1 and 2, respectively (Fig. 5). Snow both accumulated and melted during each melt period, thus the changes in SWE calculated from profiles are considered loss of SWE between surveys rather than melt alone. The maximum loss of SWE occurred in the SO for melt period 1 and in FO for melt period 2 though it is important to note that some of the SO area had a SWE of zero on May 4 and all of the snow was melted on the entire south facing slope by May 18. Because of this the minimum loss of SWE for melt period 2 occurred on the SO with a value of zero (Fig. 5(c)). The minimum loss of SWE for melt period 1 was 11 mm occurring in both FO-2 and the FC-1 plot (Fig. 5).

Comparing the survey plots to the SNOTEL station, mean SWE was generally less than the SNOTEL with the overall mean always less (Fig. 6). The SNOTEL site measured at least 120 mm greater than the overall mean for all survey dates. Daily SWE measurements at the SNOTEL station additionally captured the accumulation that occurred during both melt periods (10 mm on April 30, 5 mm on May 2, 3 mm on May 12, and 15 mm on May 13) and a maximum daily melt during the observation period of 51 mm (Fig. 7). For the three survey dates SWE measurements resulted in increasing coefficients of variation with time (Fig. 6). The FO survey plots resulted in the largest measured SWE and the lowest coefficients of variation. The highest coefficient of variation was observed in FD on May 18; these survey plot types were also the locations with the shallowest measurable snowpack of all locations. As the coefficient of variation for SWE increased through the melt season, the SNOTEL station observed losses increased in error when compared with the mean of all survey plots (Fig. 8). When observing the loss of SWE for the two melt periods, the SNOTEL station and overall mean loss of SWE differed by 7 mm (4%) for melt period 1 and 41 mm (20%) for melt period 2. The coefficients of variation for loss of SWE observed an increase for SO and FD through time and a decrease for FO and FC (Fig. 8). However, this is the result of numerous zero changes during melt period 2 when the south aspect had patches of no snow for both surveys. Generally, the overall coefficient of variation for loss of SWE remained roughly the same.

Most plot types display bimodal distributions of SWE (Figs. 9(a)‒9(d)). FD and FC histograms display a slight bimodal tendency on April 20, and normal distribution May 4 and 18. The histograms for loss of SWE, however, display a bimodal pattern for all cover types, though SO and FC plots resulted in trimodal loss of SWE distributions during melt period 2 (Fig. 9). These patterns show how variable the loss of SWE can be during spring snowmelt and that the distribution for these losses differs from the distribution of SWE.

SWE variograms show the different ranges of correlation for survey plot types in space for all three surveys (Fig. 10). Approximately half of the survey plots (6 of 11) show a decreasing range with time while roughly a third (4 of 11) show increasing range and one plot shows no change in range with time. In general, the range observed in SWE variograms show correlation lengths between two and five meters with six of the 31 observations having one to two meter correlation lengths and 11 observations having five to nine meter correlation lengths (Fig. 10). Most of the larger correlation lengths were observed in FD plots (Figs. 10(g)‒10(k)). Survey plots of similar cover type show no patterns in range increasing or decreasing with time and consistently resulted in differing range magnitudes between survey plots. Additionally, survey plots of similar aspect and cover type with surveys conducted in the same directional orientation show large differences in range (e.g., Fig. 10(h) and Fig. 10(i)) indicating that influences other than aspect and canopy cover are causing this variability.

Variograms for the loss of SWE resulted in more apparent changes with time (Fig. 11). FO survey plots tended to decrease in range with time (Figs. 11(a)‒11(c)) whereas FD tended towards increasing range with time (Figs. 11(g)‒11(k)). One plot (FD-4) survey resulted in no sill in the variogram for the melt period 1 whereas the melt period 2 resulted in a clear sill with a relatively low range of 1.2 m (Fig. 11(j)). This was the only plot that displayed this drastic of a difference between melt periods. The majority of correlation lengths observed for loss of SWE were between two and four meters with four of the 22 observations below two meters, and only three observed correlation lengths above five meters (Fig. 11). The correlation lengths for loss of SWE were lower than the correlation lengths of SWE indicating higher variability in melt relative to SWE.

Discussion

Results of this study indicate the variability in SWE and loss of SWE for different aspect and canopy covers during spring snowmelt. Survey plots in relatively close proximity display different results in both SWE and loss of SWE indicating that using a single plot to be representative of SWE or loss of SWE for a particular aspect and canopy cover type may be inaccurate, especially since SWE variability increased with time during the melt season (Fig. 6). The variability in loss of SWE was observed to be similar for the two melt periods and the SNOTEL station reasonably captured the overall mean loss of SWE for the surrounding area, though with increasing error as the melt season progressed (Fig. 8). Correlation lengths in SWE and loss of SWE, as shown by variograms (Fig. 9 and Fig. 10), resulted in lesser correlation lengths for loss of SWE relative to SWE indicating that patterns in SWE may not be appropriate for predicting the spatial patterns of melt. However, the two to four meter correlation length of loss of SWE is likely captured by the SNOTEL pillow that has a nine square meter footprint explaining the reasonable comparison with the overall mean loss observed in the study plots during melt period 1.

When comparing two study plots of similar slope aspect and canopy type, it is clear that more than these two factors will influence SWE and melt. Previous studies have shown that wind redistribution will have a strong influence on the accumulation of snow (Litaor et al., 2008; Broxton et al., 2015). For snowmelt, it is important to also consider sky view factor (SVF) of nearby terrain and/or tree canopies (López-Moreno and Latron, 2008; Musselman et al., 2012). The SVF can have a strong effect on the amount of incoming shortwave radiation at a location that will be influenced by the time of year and resulting path of the sun across the sky (Musselman et al., 2012). Vegetation canopy will further add complexity to the energy balance from emitting longwave radiation to surrounding areas (Winkler et al., 2005; Staples et al., 2006; Varhola et al., 2010). Canopy interception during the survey period needs to also be considered as a potential source of error. However, these errors introduced are minimal with only 15 mm of SWE falling during melt period 1 and 18 mm during melt period 2 (Fig. 7). With the overall mean loss of SWE for each of these periods being 181 mm and 208 mm, snowfall represents 8% and 9% of melt, respectively (Fig. 8). The largest error from canopy interception would likely occur in the FC plot where these amounts of precipitation could result in errors that are 15% and 8% of mean calculated loss of SWE for melt period 1 and melt period 2, respectively (Fig. 8). These reasons, in addition to the slope and canopy cover type, will affect the snowmelt dynamics for the plots surveyed and future investigations should observe shorter time periods without precipitation to improve understanding of snowmelt variability.

The differences among plots is further highlighted in the bimodal distributions of SWE for many of the plot types (Fig. 9), though not for all survey dates. This indicates that while SWE of similar plots may show similar patterns the loss of SWE during melt is variable enough to alter SWE distribution patterns. This is shown in the loss of SWE histograms showing at least bimodal distribution with some trimodal distribution patterns (Fig. 9). Furthermore, the correlation lengths for loss of SWE in FO plots decreased with time indicating an increase in variability whereas the correlation lengths for FD increased with time indicating a decrease in variability and homogenizing of the snowpack relative to open areas. However, the variability among plots of similar aspect and canopy type displays the complexity of snow accumulation and melt processes and that using single point/plot measurements to be representative of that aspect and land cover type may not accurately scale up for hydrological modeling or flux estimates.

At this location the SNOTEL site did capture the overall mean loss of SWE well, particularly for melt period 1 with an error of only 7 mm (Fig. 8). The SNOTEL station SWE measurements were consistently larger than most of the surveyed plots (Fig. 6) as a result of the location having no canopy for intercepting falling snow but nearby trees to provide solar shading and wind deposition (Broxton et al., 2015). The DL SNOTEL site has some solar shading during the time of this study primarily from nearby deciduous stands. This provides less shading than coniferous trees and more than open areas far from trees such as the south aspect slope. During melt period 1 this resulted in reasonably capturing the overall mean loss of SWE for the surrounding areas. It is unclear how well the SNOTEL station would perform in capturing the loss of SWE later into the spring snowmelt season as the snow becomes more “patchy”, though the error did increase to 41 mm during melt period 2 suggesting error would increase with time. For modeling purposes more seasons and melt conditions are required in order to draw stronger conclusions concerning the representative capability of this SNOTEL station to capture losses in SWE during spring. However, this study provides evidence that a SNOTEL station is capable of capturing the mean loss in SWE for an area larger than the footprint of the pillow during spring snowmelt at the weekly time scale.

Spatial patterns of spring snowmelt observed in this study show a range of correlation lengths for both SWE and loss of SWE over time. In general, loss of SWE correlation lengths were observed to be less than the five to seven meters found by other studies in alpine regions (Sommerfeld et al., 1994; Williams et al., 1999)(Fig. 10). Additionally, the correlation lengths of the loss of SWE was generally less than the correlation lengths of SWE (Fig. 11) indicating the complex nature of a melting snowpack. Aerial photographs used bySommerfeld et al. (1994) and Williams et al. (1999) were analyzed by observing darker regions that are indicative of depressions in the snowpack that meltwater flows into causing increased metamorphism and thus larger grain sizes. These darker areas will have a lower albedo, absorbing more shortwave radiation and propagating the preferential melt pattern. In subalpine terrain, however, vegetative shading will cause variability in incoming shortwave radiation in addition to the variable albedo. Furthermore, the surface roughness can have high variability (Fassnacht et al., 2009) that will result in a range of incidence angles across the surface of the snowpack with incoming shortwave radiation. These added complexities in mixed canopy subalpine terrain results in the correlation lengths observed in SWE and loss of SWE measurements. The above mentioned shortwave shading and change in incidence angles will change with time as the path of the sun across the sky moves during the spring season (Musselman et al., 2012) and so future investigations may benefit from observing shorter time scales than the two week melt periods in this study as well as over the entire melt season until all the snow has disappeared over the entire study area.

The differences in GPR wave velocity for plots can be explained by differences in dry snow density and liquid water content (Heilig et al., 2009; Mitterer et al., 2011; Schmid et al., 2015). In this study, a weaker reflection can be seen from the ground surface in some locations (Fig. 2). This is the result of melt and flow paths of liquid water either throughout the entire snowpack or accumulating near the snow-soil-interface (Granlund et al., 2009; Gusmeroli and Grosse, 2012). No field measurements were taken to observe the liquid water content of the snowpack at any location so the differences in this parameter was not assessed for spatial patterns. However, future investigations could benefit from this spatial analysis for improving upon SWE estimates for an area.

Other sources of uncertainty in these observations include the number of density and depth measurements used for calibration. However, when observing the differences in density a maximum difference of 15% is shown for plot types (Fig. 3). This magnitude of difference is slightly larger than differences observed within individual snow pits (up to 10%) (Webb and Fassnacht, 2016) suggesting that the assumptions made in this study introduces minimal errors on the order of five to ten percent that is considered reasonable. Future investigations will benefit from taking more depth measurements along each GPR transect or utilizing light detection and ranging (LiDAR). Furthermore, the use of differential GPS systems integrated with the GPR system could greatly improve accuracy of locating each measurement along transects geospatially. The combination of differential GPS and LiDAR would also allow for a more thorough measurement strategy that collects more data within each plot to further analyze controls on SWE and melt variability in subalpine environments.

Understanding the variability of snowmelt within an area can be important for hydrological modeling and water resource management. This study shows that the variability of snowmelt is different from SWE distribution because of a number of factors including canopy cover, wind shielding, slope, and aspect similar to driving parameters of variabitily during accumulation processes. However, snowmelt has the added complexity of shortwave radiation variability from the SVF and radiation incidence angles from the surface roughness. The variable melt that is different from that of SWE can partially explain variable signals in soil moisture response observed in complex subalpine terrain (Bales et al., 2011; Harpold et al., 2015; Webb et al., 2015). Future studies that observe the signal of snowmelt in soil moisture and estimate ground water recharge could benefit from applying similar methods as shown in this study to assess the variability of snowmelt through non-destructive methods.

Conclusions

The results from this investigation show the high variability for SWE and loss of SWE during spring snowmelt for varying slope, aspect and land cover types in a subalpine mixed canopy forest. In the observed study site locations with similar aspect and canopy type show different spatial patterns, even if separated by as little as 30 m. The variability in spatial patterns of SWE was shown to increase with time during the spring snowmelt season whereas the variability in loss of SWE remained similar for both melt periods observed. These observations showed a pattern of melt variability increasing in flat open areas and decreasing in flat areas with deciduous cover, though further investigations are necessary. During these melt periods for all aspect and land cover types SWE was shown to have correlation lengths between two and five meters and the loss of SWE showed correlation lengths between two and four meters. These results show lesser correlation lengths than those found in alpine regions due to the added complexity of snowmelt from the sky view factor of tree canopies. The SNOTEL station regularly measured higher SWE relative to the majority of survey plots, but did capture the loss of SWE reasonably well when compared to the overall mean of surveyed plots during melt period 1. Survey plots of similar aspect and land cover type varied in SWE and loss of SWE, indicating that care should be taken if measuring a single location and scaling up for hydrological modeling purposes. Hydrological investigations can benefit from the non-destructive nature of GPR methods that will improve SWE surveys and melt flux estimates for regions in the future. Recommendations for future investigations include utilizing differential GPS and LiDAR systems that would enable more accurate measurements at greater spatial resolution. Additionally, more snow pit observations with measured liquid water contents would enhance observations of meltwater distribution within the snowpack. These recommendations for future research will improve our understanding of the controls on melt distribution and variability in subalpine terrain that is driven by complex physical processes.

References

[1] Adam J C, Hamlet  A F, Lettenmaier  D P (2009). Implications of global climate change for snowmelt hydrology in the twenty-first century. Hydrol Processes, 23(7): 962–972

[2] Andreadis K M ,  Storck P ,  Lettenmaier D P  (2009). Modeling snow accumulation and ablation processes in forested environments. Water Resour Res, 45(5) 

[3] Bales R C, Hopmans  J W, O’Geen  A T, Meadows  M, Hartsough P C ,  Kirchner P ,  Hunsaker C T ,  Beaudette D  (2011). Soil moisture response to snowmelt and rainfall in a Sierra Nevada mixed-conifer forest. Vadose Zone J, 10(3): 786–799

[4] Bales R C, Molotch  N P, Painter  T H, Dettinger  M D, Rice  R, Dozier J  (2006). Mountain hydrology of the western United States. Water Resour Res, 42(8): W08432

[5] Blöschl G (1999). Scaling issues in snow hydrology. Hydrol Processes, 13(14–15): 2149–2175

[6] Blöschl G, Kirnbauer  R (1992). An analysis of snow cover patterns in a small alpine catchment. Hydrol Processes, 6(1): 99–109

[7] Broxton P D, Harpold  A A, Biederman  J A, Troch  P A, Molotch  N P, Brooks  P D (2015). Quantifying the effects of vegetation structure on snow accumulation and ablation in mixed-conifer forests. Ecohydrology, 8(6): 1073–1094

[8] Cao J, Liu  C, Zhang W  (2012). Response of rock-fissure seepage to snowmelt in Mount Taihang slope-catchment, North China. Water Sci Technol, 67(1): 124–130

[9] Clark M P, Hendrikx  J, Slater A G ,  Kavetski D ,  Anderson B ,  Cullen N ,  Kerr T, Örn Hreinsson  E, Woods R A  (2011). Representing spatial variability of snow water equivalent in hydrologic and land-surface models: a review. Water Resour Res, 47(7): W07539

[10] Clilverd H M, White  D M, Tidwell  A C, Rawlins  M A (2011). Sensitivity of northern groundwater recharge to climate change: a case study in Northwest Alaska. Journal of the American Water Resources Association, 47(6): 1228–1240 

[11] Cline D, Yueh  S, Chapman B ,  Stankov B ,  Gasiewski A ,  Masters D ,  Elder K ,  Kelly R ,  Painter T H ,  Miller S ,  Katzberg S ,  Mahrt L  (2009). NASA Cold Land Processes Experiment (CLPX 2002/03): airborne remote sensing. J Hydrometeorol, 10(1): 338–346

[12] Clow D W (2010). Changes in the timing of snowmelt and streamflow in Colorado: a response to recent warming. J Clim, 23(9): 2293–2306

[13] Daly S F, Davis  R, Ochs E ,  Pangburn T  (2000). An approach to spatially distributed snow modelling of the Sacramento and San Joaquin basins, California. Hydrol Processes, 14(18): 3257–3271

[14] Ebel B A, Hinckley  E S, Martin  D A (2012). Soil-water dynamics and unsaturated storage during snowmelt following wildfire. Hydrol Earth Syst Sci, 16(5): 1401–1417

[15] Eiriksson D, Whitson  M, Luce C H ,  Marshall H P ,  Bradford J ,  Benner S G ,  Black T ,  Hetrick H ,  McNamara J P  (2013). An evaluation of the hydrologic relevance of lateral flow in snow at hillslope and catchment scales. Hydrol Processes, 27(5): 640–654

[16] Elder K, Cline  D, Liston G E ,  Armstrong R  (2009). NASA Cold Land Processes Experiment (CLPX 2002/03): field measurements of snowpack properties and soil moisture. J Hydrometeorol, 10(1): 320–329

[17] Elder K, Dozier  J, Michaelsen J  (1991). Snow accumulation and distribution in an Alpine watershed. Water Resour Res, 27(7): 1541–1552

[18] Fang S, Xu  L, Zhu Y ,  Liu Y, Liu  Z, Pei H ,  Yan J, Zhang  H (2015). An integrated information system for snowmelt flood early-warning based on internet of things. Inf Syst Front, 17(2): 321–335

[19] Fassnacht S R ,  Cherry M L ,  Venable N B H ,  Saavedra F  (2016). Snow and albedo climate change impacts across the United States Northern Great Plains. Cryosphere, 10(1): 329–339

[20] Fassnacht S R ,  Derry J E  (2010). Defining similar regions of snow in the Colorado River Basin using self-organizing maps. Water Resour Res, 46(4): W04507

[21] Fassnacht S R ,  Hultstrand M  (2015). Snowpack variability and trends at long-term stations in northern Colorado, USA. International Association of Hydrological Sciences, 92: 1–6 

[22] Fassnacht S R ,  Williams S R ,  Corrao M V  (2009). Changes in the surface roughness of snow from millimetre to metre scales. Ecol Complex, 6(3): 221–229

[23] Flint A L, Flint  L E, Dettinger  M D (2008). Modeling soil moisture processes and recharge under a melting snowpack. Vadose Zone J, 7(1): 350–357

[24] Granlund N, Lundberg  A, Feiccabrino J ,  Gustafsson D  (2009). Laboratory test of snow wetness influence on electrical conductivity measured with ground penetrating radar. Hydrol Res, 40(1): 33–44

[25] Graybeal D, Leathers  D (2006). Snowmelt-related flood risk in Appalachia: first estimates from a historical snow climatology. J Appl Meteorol Climatol, 45(1): 178–193

[26] Gusmeroli A, Grosse  G (2012). Ground penetrating radar detection of subsnow slush on ice-covered lakes in interior Alaska. Cryosphere, 6(6): 1435–1443

[27] Harpold A, Brooks  P, Rajagopal S ,  Heidbuchel I ,  Jardine A ,  Stielstra C  (2012). Changes in snowpack accumulation and ablation in the intermountain west. Water Resour Res, 48(11): W11501

[28] Harpold A A, Biederman  J A, Condon  K, Merino M ,  Korgaonkar Y ,  Nan T, Sloat  L L, Ross  M, Brooks P D  (2014). Changes in snow accumulation and ablation following the Las Conchas Forest Fire, New Mexico, USA. Ecohydrology, 7(2): 440–452

[29] Harpold A A, Molotch  N P, Musselman  K N, Bales  R C, Kirchner  P B, Litvak  M, Brooks P D  (2015). Soil moisture response to snowmelt timing in mixed-conifer subalpine forests. Hydrol Processes, 29(12): 2782–2798

[30] Heilig A, Mitterer  C, Schmid L ,  Wever N ,  Schweizer J ,  Marshall H P ,  Eisen O  (2015). Seasonal and diurnal cycles of liquid water in snow-Measurements and modeling. J Geophys Res Earth Surf, 120(10): 2139–2154

[31] Heilig A, Schneebeli  M, Eisen O  (2009). Upward-looking ground-penetrating radar for monitoring snowpack stratigraphy. Cold Reg Sci Technol, 59(2‒3): 152–162

[32] Jencso K, McGlynn  B, Gooseff M ,  Wondzell S ,  Bencala K ,  Marshall L  (2009). Hydrologic connectivity between landscapes and streams: transferring reach-and plot-scale understanding to the catchment scale. Water Resour Res, 45(4): W04428

[33] Jencso K G, McGlynn  B L (2011). Hierarchical controls on runoff generation: topographically driven hydrologic connectivity, geology, and vegetation. Water Resour Res, 47(11): W11527

[34] Johnson J B, Schaefer  G L (2002). The influence of thermal, hydrologic, and snow deformation mechanisms on snow water equivalent pressure sensor accuracy. Hydrol Processes, 16(18): 3529–3542

[35] Litaor M I, Williams  M, Seastedt T R  (2008). Topographic controls on snow distribution, soil moisture, and species diversity of herbaceous alpine vegetation, Niwot Ridge, Colorado. J Geophys Res, 113(G2): G02008

[36] Liu F, Williams  M W, Caine  N (2004). Source waters and flow paths in an alpine catchment, Colorado Front Range, United States. Water Resour Res, 40(9): W09401

[37] López-Moreno J ,  Fassnacht S ,  Heath J ,  Musselman K ,  Revuelto J ,  Latron J ,  Moran-Tejeda E ,  Jonas T  (2013). Small scale spatial variability of snow density and depth over complex alpine terrain: implications for estimating snow water equivalent. Adv Water Resour, 55: 40–52

[38] López-Moreno J I ,  Fassnacht S R ,  Beguería S ,  Latron J B P  (2011). Variability of snow depth at the plot scale: implications for mean depth estimation and sampling strategies. Cryosphere, 5(3): 617–629

[39] López-Moreno J I ,  Latron J  (2008). Influence of canopy density on snow distribution in a temperate mountain range. Hydrol Processes, 22(1): 117–126

[40] Lundberg A, Ala-Aho  P, Eklo O ,  Klove B ,  Kvaerner J ,  Stumpp C  (2016). Snow and frost: implications for spatiotemporal infiltration patterns–A review. Hydrol Processes, 30(8): 1230–1250

[41] Lundberg A, Richardson-Naslund  C, Andersson C  (2006). Snow density variations: consequences for ground-penetrating radar. Hydrol Processes, 20(7): 1483–1495

[42] Lundquist J, Cayan  D (2002). Seasonal and spatial patterns in diurnal cycles in streamflow in the western United States. J Hydrometeorol, 3(5): 591–603

[43] Lundquist J, Dettinger  M, Cayan D  (2005). Snow-fed streamflow timing at different basin scales: case study of the Tuolumne River above Hetch Hetchy, Yosemite, California. Water Resour Res, 41(7): W07005

[44] Magnusson J, Kobierska  F, Huxol S ,  Hayashi M ,  Jonas T ,  Kirchner J W  (2014). Melt water driven stream and groundwater stage fluctuations on a glacier forefield (Dammagletscher, Switzerland). Hydrol Processes, 28(3): 823–836

[45] McNamara J P, Chandler  D, Seyfried M ,  Achet S  (2005). Soil moisture states, lateral flow, and streamflow generation in a semi-arid, snowmelt-driven catchment. Hydrol Processes, 19(20): 4023–4038

[46] Mitterer C, Heilig  A, Schweizer J ,  Eisen O  (2011). Upward-looking ground-penetrating radar for measuring wet-snow properties. Cold Reg Sci Technol, 69(2‒3): 129–138

[47] Moeser D, Mazzotti  G, Helbig N ,  Jonas T  (2016). Representing spatial variability of forest snow: implementation of a new interception model. Water Resour Res, 52(2): 1208–1226

[48] Molotch N P, Brooks  P D, Burns  S P, Litvak  M, Monson R K ,  McConnell J R ,  Musselman K  (2009). Ecohydrological controls on snowmelt partitioning in mixed-conifer sub-alpine forests. Ecohydrology, 2(2): 129–142

[49] Molotch N P, Colee  M T, Bales  R C, Dozier  J (2005). Estimating the spatial distribution of snow water equivalent in an alpine basin using binary regression tree models: the impact of digital elevation data and independent variable selection. Hydrol Processes, 19(7): 1459–1479

[50] Molotch N P, Meromy  L (2014). Physiographic and climatic controls on snow cover persistence in the Sierra Nevada Mountains. Hydrol Processes, 28(16): 4573–4586

[51] Musselman K N ,  Molotch N P ,  Brooks P D  (2008). Effects of vegetation on snow accumulation and ablation in a mid-latitude sub-alpine forest. Hydrol Processes, 22(15): 2767–2776

[52] Musselman K N ,  Molotch N P ,  Margulis S A ,  Kirchner P B ,  Bales R C  (2012). Influence of canopy structure and direct beam solar irradiance on snowmelt rates in a mixed conifer forest. Agric Meteorol, 161: 46–56

[53] Mutzner R, Weijs  S, Tarolli P ,  Calaf M ,  Oldroyd H ,  Parlange M  (2015). Controls on the diurnal streamflow cycles in two subbasins of an alpine headwater catchment. Water Resour Res, 51(5): 3403–3418

[54] Previati M, Godio  A, Ferraris S  (2011). Validation of spatial variability of snowpack thickness and density obtained with GPR and TDR methods. J Appl Geophys, 75(2): 284–293

[55] Rice R, Bales  R C, Painter  T H, Dozier  J (2011). Snow water equivalent along elevation gradients in the Merced and Tuolumne River basins of the Sierra Nevada. Water Resour Res, 47(8): W08515

[56] Richer E E, Kampf  S K, Fassnacht  S R, Moore  C C (2013). Spatiotemporal index for analyzing controls on snow climatology: application in the Colorado Front Range. Phys Geogr, 34(2): 85–107

[57] Schmid L, Koch  F, Heilig A ,  Prasch M ,  Eisen O ,  Mauser W ,  Schweizer J  (2015). A novel sensor combination (upGPR-GPS) to continuously and nondestructively derive snow cover properties. Geophys Res Lett, 42(9): 3397–3405

[58] Sexstone G A, Fassnacht  S R (2014). What drives basin scale spatial variability of snowpack properties in northern Colorado? Cryosphere, 8(2): 329–344

[59] Seyfried M S, Grant  L E, Marks  D, Winstral A ,  McNamara J  (2009). Simulated soil water storage effects on streamflow generation in a mountainous snowmelt environment, Idaho, USA. Hydrol Processes, 23(6): 858–873

[60] Singh K K, Datt  P, Sharma V ,  Ganju A ,  Mishra V D ,  Parashar A ,  Chauhan R  (2011). Snow depth and layer interface estimation using Ground Penetrating Radar. Curr Sci, 100(10): 1532–1539

[61] Sommerfeld R A ,  Bales R C ,  Mast A (1994). Spatial statistics of snowmelt flow-data from lysimeters and aerial Photos. Geophys Res Lett, 21(25): 2821–2824

[62] Staples J M, Adams  E E, Slaughter  A E, McKittrick  L R (2006). Slope scale modeling of snow surface temperature in topographically complex terrain. In: Proceedings of 2006 International Snow Science Workshop: 806–814

[63] Storck P, Lettenmaier  D P, Bolton  S M (2002). Measurement of snow interception and canopy effects on snow accumulation and melt in a mountainous maritime climate, Oregon, United States. Water Resour Res, 38(11): 1223

[64] USGS (2015). 3DEP products and services: The national Map, 3D Elevation Program Web page. accessed (November, 2015), edited, 

[65] Varhola A, Coops  N C, Weiler  M, Moore R D  (2010). Forest canopy effects on snow accumulation and ablation: an integrative review of empirical results. J Hydrol (Amst), 392(3-4): 219–233

[66] Webb R W, Fassnacht  S R (2016). Snow density, snow depth, and soil moisture at Dry Lake study site in northern Colorado in 2014. Colorado State University, 

[67] Webb R W, Fassnacht  S R, Gooseff  M N (2015). Wetting and drying variability of the shallow subsurface beneath a snowpack in California’s Southern Sierra Nevada. Vadose Zone J, 14(8): doi: 10.2136/vzj2014.12.0182

[68] Williams C J, McNamara  J P, Chandler  D G (2009a). Controls on the temporal and spatial variability of soil moisture in a mountainous landscape: the signature of snow and complex terrain. Hydrol Earth Syst Sci, 13(7): 1325–1336

[69] Williams M W, Seibold  C, Chowanski K  (2009b). Storage and release of solutes from a subalpine seasonal snowpack: soil and stream water response, Niwot Ridge, Colorado. Biogeochemistry, 95(1): 77–94

[70] Williams M W, Sommerfeld  R, Massman S ,  Rikkers M  (1999). Correlation lengths of meltwater flow through ripe snowpacks, Colorado Front Range, USA. Hydrol Processes, 13(12–13): 1807–1826

[71] Winkler R D, Spittlehouse  D L, Golding  D L (2005). Measured differences in snow accumulation and melt among clearcut, juvenile, and mature forests in southern British Columbia. Hydrol Processes, 19(1): 51–62

[72] Zhao Q, Liu  Z, Ye B ,  Qin Y, Wei  Z, Fang S  (2009). A snowmelt runoff forecasting model coupling WRF and DHSVM. Hydrol Earth Syst Sci, 13(10): 1897–1906

RIGHTS & PERMISSIONS

Higher Education Press and Springer-Verlag Berlin Heidelberg