Meltwater Storage in the firn of Kaskawulsh Glacier, Yukon Territory, Canada

In recent years, the analysis of firn in Greenland, Svalbard, and other high Arctic regions has contributed to the 10 understanding of meltwater retention in firn and its importance to measurements of glacier mass balance. This has provided insight into firn densification processes and meltwater retention. Changes in these attributes can also provide insight into meteorological variability and climate trends. In spring 2018, two firn cores (21 m and 36 m in length) were extracted from the accumulation zone of Kaskawulsh Glacier, St. Elias Mountains, Yukon. The cores were analyzed for ice layer stratigraphy, density, and glaciochemical time series (oxygen isotopes and major ions). Meltwater percolation and refreezing 15 events were evident in the cores. The quantity of ice layers, the presence of liquid water at 34.5 m depth, interpreted as a perennial firn aquifer (PFA), and the altered isotopic and glaciochemical signature all indicate this process. This melt resulted in an estimated surface lowering of 10 ± 0.8 cm/yr between 2005 and 2018. The information gleaned from Kaskawulsh Glacier supports the need for improved and field-validated density assumptions for geodetic mass balance methods. 20


Introduction
With the increasing affects of climate change and the need for understanding glacier and ice sheet melt rates, geodetic methods are useful for indirect measurements of mass balance. Based on repeat altimetry, geodetic approaches to mass balance monitoring rely on several assumptions. Estimates must be made of the density of snow, firn, and ice at the sampling 25 location, with the additional assumption that these densities remain unchanged between the two measurement dates. With large changes in elevation that can occur over decadal timescales, changes in density of the firn may have little impact on the mass balance calculations from geodetic measurements (Moholdt et al., 2010a). However, over shorter time periods (years), this may not be true (Moholdt et al., 2010b). Meltwater percolation and refreezing can significantly change the firn density profile and density of the accumulation zone of a glacier . These modifications of firn density can 30 introduce large uncertainties when using geodetic techniques to determine glacier mass balance. For example, Moholdt et al. (2010a) determined mass balance on a Svalbard glacier to be −4.3±1.4 Gt/yr, with the large uncertainty attributed to limited https://doi.org/10.5194/tc-2020-119 Preprint. Discussion started: 25 May 2020 c Author(s) 2020. CC BY 4.0 License. knowledge of the snow and firn density and their spatial and temporal variation. By altering the density and causing surface lowering, meltwater percolation, refreezing, and liquid water storage all complicate the interpretation of geodetic mass balance data. 35 Warming firn can result in increased meltwater production and altered firn densification processes. Initially, melt can round the snow grains and increase the snowpack density. If the surface continues to melt, the meltwater can percolate into the firn and refreeze as ice layers or lenses. On temperate glaciers, meltwater that percolates below the winter cold layer often will not refreeze, and may thus form a firn aquifer (Munneke et al, 2014). These internal accumulation processes can 40 significantly affect the density of the firn. Once ice layers or firn aquifers form, this also effects how meltwater percolates through the firn thereafter. Due to the spatial heterogeneity of meltwater retention, percolation, and refreezing processes, there are still many gaps in knowledge on how to model these processes and subsequently estimate firn density in areas where these processes occur (van As et al., 2016).

45
Meltwater retention in firn is also important for estimating glacial runoff contributions to sea level rise. Recently, there have been extensive studies investigating the meltwater and refreezing occurring in southern Greenland (Humphrey et al., 2012;Harper et al., 2012;De La Peña et al., 2015;Macferrin et al., 2019), Arctic Canada Gascon et al., 2013), and Svalbard (Christianson et al., 2015). Using ground and airborne radar as well as ice cores, researchers have estimated the volume of water storage in the percolation zone of the Greenland Ice Sheet and have investigated the 50 interactions of the meltwater and ice and their subsequent effects on hydrology and run-off mechanisms (e.g., Harper et al., 2012;Koenig et al., 2014;Machguth et al., 2016). Increased melt rates in the Greenland percolation zone may lead to the expansion of low-permeability ice layers, causing run-off in Greenland to increase and expediting the movement of water from the ice sheet to the ocean (Macferrin et al., 2019). These processes are not limited to the Arctic; it is necessary to investigate these properties on mountain glaciers as well, for improved estimates and models of glacier mass balance and 55 associated sea-level rise.
In this study two firn cores were retrieved on Kaskawulsh Glacier, St. Elias Icefields, Yukon, Canada. The firn cores were analyzed for density and the effects of meltwater percolation and refreezing. This included physical descriptions and quantification of ice layers and lenses, as well as analyses of the stable isotopes δ 18 O and δD and major ions. Foy et al., 60 (2011) recorded density up to 8.1m depth in 2007 in the accumulation zone of Kaskawulsh Glacier. However, recent and detailed documentation of firn density in the St. Elias Icefields is lacking. The region was also the location of numerous pioneering glaciological studies in the 1960s, documented in a series of Icefield Ranges Research Reports (e.g., IRRP 1963, 1965, 1972. These early studies provide an estimate of firn density within the upper 15 m from a core close to our drilling site, providing an estimate of changes over the past ~50 years. 65 https://doi.org/10.5194/tc-2020-119 Preprint. Discussion started: 25 May 2020 c Author(s) 2020. CC BY 4.0 License.

Study Area
We retrieved two firn cores in May 2018 from the upper accumulation zone of Kaskawulsh Glacier in the St. Elias Mountains (Figure 1), near the icefield divide. The St. Elias Mountains are located in the southwest corner of Yukon Territory, Canada (Figure 1). The peaks in this mountain range are predominantly between 300 -2200 m above sea level 70 (a.s.l.), with several notable exceptions of mountains higher than 3200 m, including Mount Logan at 5959 m a.s.l. (Foy et al., 2011). The St. Elias range is home to the largest icefield outside of the polar regions, with an area of 46,000 km 2 (Berthier et al., 2010). Like most high northern latitudes, Yukon Territory is warming (Foy et al., 2011;Williamson et al., 2020) and there is a need to better understand and document the impacts of climate change on the St. Elias Icefield.

75
Kaskawulsh Glacier is a large temperate valley glacier located on the eastern side of the St. Elias Mountains within the Donjek Range, and is approximately 70 km long and 3 -4 km wide. The drill site is located on the upper north arm of the glacier in the accumulation zone (60.78°N, 139.63°W), at an elevation of ~2640 m a.s.l. Previous studies at Kaskawulsh Glacier include an analysis of volume change over time using satellite imagery (Foy et al., 2011). There were also several reports in the early 1960s that documented various glaciological characteristics and processes occurring in the St. Elias 80 Mountains (Grew and Mellor, 1966;Wood, 1963). A weather station (Copland Weather Station, Figure 1) is 12 km southwest of the drill site (60.70°N, 139.80°W, 2600 m a.s.l.) has been in service since 2013 (Williamson et al., 2020). The temperature sensor is an Onset S-THB-M002 with an accuracy of ± 0.2° C, mounted inside an Onset RS3 naturally vented solar radiation shield. There was a recent expedition (June 2018) at a nearby location, Kreutz Drill Site ( Figure 1) (Kreutz, personal comm.) Other studies in the region include analysis of ice cores from the Eclipse Icefield, located 12 km northwest 85 of our drill site (e.g., Zdanowicz et al., 2014 andYalcin et al., 2006).

Field Methods
Two 8-cm diameter cores were drilled using an ECLIPSE ice core drill (Icefield Instruments, Yukon Territory, Canada). The seasonal snow had a depth of 4.4 m. With a starting depth of 2 m below the snow surface, Core 1 was 34.6 m long and reached a depth of 36.6 m, and Core 2 was 19.5 m long and reached a depth of 21.5 m. The two cores were drilled 60 cm 90 apart, and core stratigraphy and density were recorded in the field. At a depth of 34.5 m during the drilling of Core 1, liquid water became evident; drilling was stopped at a depth of 36.6 m to avoid the risk of the drill freezing in the hole. A complementary 1 m deep snow pit was dug several meters away, allowing for detailed surface density measurements and sample collection for chemical and isotopic analysis.

95
Once the cores were retrieved, the presence of ice layers, lenses and melt-affected firn was logged and the stratigraphic character, depth, and thickness were recorded. When an ice horizon extended across the entire width of the core, it was https://doi.org/10.5194/tc-2020-119 Preprint. Discussion started: 25 May 2020 c Author(s) 2020. CC BY 4.0 License. labeled as an ice layer. If the ice horizon was of more limited lateral extent, it was labeled an ice lens. Ice lenses were occasionally wedge shaped. Other types of melt-affected firn were identified by the lack of grain boundaries and the presence of air bubbles. Other noticeable features were also documented (e.g., the texture and visual appearance of the firn, 100 such as opacity, which may have indicated melt-affected firn). Each core section was then sawed into 10-cm long sections.
Each 10-cm section had its diameter measured at each end and was bagged, weighed, and assessed for the quality of the core sample and its cylindrical completeness (f). The average diameter was used to determine the volume of the core section (V).
Together with the mass of the core section, m, density was calculated: Where ρ is the density of the firn, D is the average core section diameter, L is the length, and f is the subjectively assessed fraction of completeness of the core section. For instance, if visual inspection indicated that about 5% of the core was missing (e.g., due to missing ice chips from the core dogs of the drill head), then f would be 0.95. Outliers were removed if 110 they were not physically possible (e.g., values >917 kg/m 3 or <300 kg/m 3 at depths greater than 4 m). The ice cores were transported back to Kluane Lake Research Station frozen, where the individual sub-samples were melted at ambient temperature in sample bags. Samples were then transferred into 20-ml high-density polyethylene (HDPE) vials; no headspace was allowed in the HDPE vials, and HDPE caps were sealed with parafilm. The 20-ml samples were transported in a cooler (~4C) to Alberta for chemical analysis. During transport and melting, some Core 2 sample bags were damaged, 115 such that the density could not be calculated on all Core 2 samples.

Density Analysis
In order to calculate the uncertainty of the density calculations, numerous sources of error were assessed using standard error analysis. Random and human error were taken into account and propagated through the equations. The uncertainty in f is 120 difficult to assign, but is the largest source of uncertainty in the density calculations. Three different people performed this evaluation with occasional inter-comparison of values and therefore there is subjectivity in the f values. To be conservative, the following uncertainties were assigned for f: The uncertainty of a higher f value is lower because it is more obvious when a core is of good quality; less complete cylinders are more difficult to assess, hence the greater uncertainty for values of f ≤ 0.8. Here, an upper density bound of https://doi.org/10.5194/tc-2020-119 Preprint. Discussion started: 25 May 2020 c Author(s) 2020. CC BY 4.0 License. 917 kg/m 3 (the density of pure ice) was assigned if the calculated density was larger than physically possible. The f values 130 for Core 2 were not recorded in the field. The average core quality value of Core 1 was used for Core 2, given the likelihood that the overall quality of the cores would be similar due to their very close proximity.
Ice content (i.e. ice fraction), F i , was calculated for each 10-cm section of the firn core. Here ice was defined based on its lack of air bubbles and crystalline structure, as compared to granular structure. Ice content was calculated instead of melt 135 percent, as melt percent generally assumes the meltwater remains within the net annual accumulation layer, which cannot be assumed here. The thickness of individual ice layers was summed within each 10-cm core section. On average, the ice lenses were observed to occupy 50% of the core samples; therefore their thickness was divided by two before being summed. For each core section, total ice content was divided by the length of the section, L, to give F i .

140
To understand the firn densification process in the absence of refrozen meltwater, the 'background' firn density is of interest.
For each sample, we estimate this by subtracting the mass and volume of the ice to give the firn density in the absence of ice content. We use a 30-cm moving average of ice content and density in order to smooth out a possible deviation of ±10 cm in assigning the location of the ice features within the stratigraphy. Here, ice layers and lenses were assumed to have a density of 874 kg/m 3 , based on the average density of firn-core sections that were 100% ice in Greenland (873 kg/m 3 ) and Devon Ice 145 Cap (875 kg/m 3 ) (Machguth et al., 2016;Bezeau et al., 2013). This is different than the 917 kg/m 3 upper bound used in the outlier analysis because that is the theoretical limit for pure ice, whereas this density has been determined using measured data from two other sites as well and incorporates the fact that the observed ice layers and lenses had tiny bubbles in them, indicating they were not pure ice. Each sample has a measured bulk density, ρ b , which we assume to result from a binary mixture of ice and firn, with densities ρ f and ρ i . Ice and firn fractions, F i and F f , were defined with F i + F f = 1. The 150 background firn density is then calculated following: In some cases there is no ice content (F i = 0), so ! = ! . Surface lowering associated with refreezing was calculated for each core section using the background firn density, ! , and length of the section, L. For each core sample, the 'thinning' or surface lowering associated with the ice content, ΔL, can be estimated by reverting the ice to the density of the background 155 firn, as follows: This is summed for the whole core to give the full amount of surface lowering associated with meltwater percolation and refreezing.

Stable Isotope and Ion Analysis 160
At the University of Calgary Isotope Science Laboratory (ISL), the liquid samples were pipetted into 2 ml glass isotope vials

Density 170
The first 4.4 m of both cores was observed to be snow. At a depth of 4.4 m there was significant ice crust, representing the 2017 summer melt surface. This surface snow layer was also dry. The 4.4 m of snow had an average density of 450 ± 21 kg/m 3 .
The overall pattern of density resembles a typical densification curve (Cuffey and Paterson, 2010). Even though there is 175 large scatter in the density values and a large uncertainty associated with each point, the firn density versus depth shows a logarithmic increase, as to be expected given Sorge's Law of densification (Bader, 1954) (Figure 2). To highlight this, a logarithmic curve was fitted to the density data. The density rapidly increases in the first 5 m below the snow surface.
Between 5 and 10 m depth the density varies widely, from 300 to 900 kg/m 3 . Below 10 m, firn density increases logarithmically from approximately 550 to 750 kg/m 3 . This corresponds with the second phase of densification in Sorge's 180 Law. Between 20 and 32 m depth there is a wide range of density values, from approximately 650 to 908 kg/m 3 . As a whole, Core 2 follows a similar pattern to Core 1, with an abrupt initial increase and then scattered values at greater depths.
For the specific calculations of firn density, we omit the upper 4.4 m of the core, which represented the seasonal snow accumulation from winter 2017-2018. Below this, in the upper 10 m of firn, Cores 1 and 2 have average densities of 558 ± 4 185 kg/m 3 and 571 ± 3 kg/m 3 respectively, giving an overall average density of 564 ± 2 kg/m 3 . Using both cores to 21 m depth, the upper 17 m of Kaskawulsh firn has an average density of 608 ± 2 kg/m 3 . Core 1 extends to a depth of 36.6 m, representing ~32 m of firn with an average density of 670 ± 2 kg/m 3 . Based on the high densities and liquid water content in the bottom sections of the core, we believe 670 ± 2 kg/m 3 provides an estimate of the average density of the upper 32.2 m of firn. 190 https://doi.org/10.5194/tc-2020-119 Preprint. Discussion started: 25 May 2020 c Author(s) 2020. CC BY 4.0 License.
Total ice content in Core 1 was 1.53 m or 1.34 m w.e. This is equivalent to 4.2% by volume (Table 1) and 5.8% by mass.
Using Eq. 2 to calculate the ice content over the full depth range of Core 1, we estimate a background firn density of 660 ± 6 kg/m 3 , 4.3% less than the bulk density of the firn (Table 1). Core 2 had a similar average firn density and ice content to the upper 21 m of Core 1, representing 3.8% by volume. 195 The total water equivalent (w.e.) of Core 1 is calculated to be 23.22 m w.e. Using the average accumulation rate of the Icefield Ranges Research Reports, 1.76 m w.e./yr (Wagner, 1963;Holdsworth, 1965), Core 1 represents approximately 13 years of net accumulation, from 2005 to 2018.

Stratigraphy 200
The stratigraphy of Kaskawulsh Glacier firn cores indicates numerous ice layers and melt-affected firn. The first 6 m of the stratigraphy display several small ice layers (<2.5 cm thick), interpreted as wind crusts in the snow (Figure 3). Several thick (>10 cm) ice layers are found between 6 and 26 m depth. The largest ice layer in Core 1 was 22 cm thick, found at 14.1 m.
At 26.4 m the ice layers and lenses disappeared. Below this the firn is almost entirely meltwater-affected, based on the appearance and texture of the firn, but is without the quantity of ice lens or layer content that was present in the first 25 m. 205 We interpret this section of the core as infiltration ice: water-saturated firn that has experienced refreezing. At 30 m depth the meltwater effects were absent and there were two small ice layers and an ice lens. At 30.6 m the firn was melt-affected again. From 34.5 to 36.6 m the core sections retrieved expelled liquid water as the sections were extracted from the core barrel.

210
In Core 2 there were numerous ice layers starting at a depth of 3.8 m, and below 4.4 m the firn was meltwater-affected.
There was a thick ice layer at 6.6 m that was 30 cm lower than a similar ice layer in Core 1 (at 6.3 m). There were numerous melt-affected layers between ice lenses much closer to the surface in Core 2. In Core 1 there were several ice layers at ~10 m depth, and these layers were not present in Core 2. At 14.4 m another section of the firn had numerous ice layers. This was also seen in Core 1 but is 20-30 cm deeper here. At 14.6 m the thickest ice layer was encountered (12 cm), corresponding 215 well with the thickest layer in Core 1. Between 16 and 21.5 m the firn was melt-affected. The differences between Core 1 and Core 2 stratigraphy can be attributed to horizontal variability in meltwater infiltration causing heterogeneous spacing of ice layers and lenses. This is known to occur at length scales less than 1-m (Parry et al., 2007;Harper et al., 2011).

Stable Isotopes
The range of δ 18 O values in the winter 2017-2018 snowpack was −19.5 to −29‰. Since these samples were acquired before 220 the onset of summer melt, this is assumed to reflect depositional values for fall through spring precipitation at this site, with no evidence for post-depositional modification. Below the seasonal snowpack, stable isotope variations in the firn are largely https://doi.org/10.5194/tc-2020-119 Preprint. Discussion started: 25 May 2020 c Author(s) 2020. CC BY 4.0 License. washed out. No summer peaks are preserved (i.e., values greater than -22‰), and most of the winter peaks are also washed out. Several negative peaks are partially preserved, indicated by red arrows in Figure 4. In Core 1 (Figure 4) the average δ 18 O for the 36.6-m core is −24.6 ± 0.7‰. There are several periods of homogenization where the variability is less than ~2‰. 225 Between 26.9 and 31.9 m the range is −24.8 to −23.2‰ with an average of −23.9‰, but below this the δ 18 O values begin to decrease linearly. The bottom part of the core has a range of −25.9 to −23.3‰ and an average of −24.8‰. Please note that all of the major ion data were so severely washed out that the time series are not presented here; these data were not useful for analyses, but also testify to complete flushing of the firn to 36-m depth.

Meltwater percolation and refreezing effects
Kaskawulsh Glacier has indications of widespread meltwater percolation and refreezing. This is evident in the washout of the isotopic signal. Meltwater is stored within the firn as ice, as indicated by the presence of ice layers and infiltration ice. It is also being stored in liquid form as a perennial firn aquifer (PFA). The water at depth is at least one year old, and likely the accumulation of numerous melt seasons; drilling occurred in May, prior to the onset of spring melt and any potential summer 235 replenishment of the aquifer. The solid and liquid phase storage mechanisms in the Kaskawulsh Glacier firn layer have different implications for the mass balance of the glacier. Liquid water is commonly found in the temperate firn of low-and mid-latitude mountain glaciers and has played an important role in meltwater storage and glacier hydrology (Fountain and Walder, 1989;Scheider, 1999). Depending on the melt, firn aquifer thickness, and temperature of the firn, the storage of liquid water at the firn-ice interface delays runoff from hours to weeks or longer (Jansson et al., 2003) and can account for as 240 much as 64% of internal accumulation, as found in Alaska and Sweden (Schneider, 1999;Trabant and Mayo, 1985).
The water in the firn aquifer on Kaskawulsh Glacier is likely to be flowing, redistributing mass. The drill site was located high in the glacier's accumulation zone, with a gently sloping surface (<0.6˚) resulting in a subtle hydraulic gradient. The low concentration of mobile ions in the saturated firn also implies that ions are being transported away from the site. If they 245 were remaining within the firn column, there would likely be greater ion concentrations near the bottom of the core. We likely drilled into the top of the water table of the firn aquifer. There may be flow along the firn-ice interface down slope, as well as possible Darcian flow within the firn aquifer itself. It may also enter the englacial hydrological system and eventually reach the bed of the glacier; for example research on Storglaciären showed that there might be a subglacial aquifer fed from the firn area (Jansson et al., 2003). 250 The liquid-phase meltwater retention on Kaskawulsh is similar to the firn aquifers found in the high-accumulation areas of southern Greenland and Svalbard (e.g., Miège et al., 2016;Christianson et al., 2015) and different than the water-saturated https://doi.org/10.5194/tc-2020-119 Preprint. Discussion started: 25 May 2020 c Author(s) 2020. CC BY 4.0 License. layers commonly found on temperate glaciers. Firn aquifers that have been studied on temperate mountain glaciers typically have a saturated layer close to the surface (for example, 5 m below the surface at Storglaciären), have active discharge and 255 recharge processes (Fountain and Walder, 1998;Schneider, 1999), and appear to experience seasonal drainage over the winter months (Fountain, 1989(Fountain, , 1996Jansson et al., 2003), likely due to high slope gradients. In 2012, "water-saturated" firn was found at 40-m depth in an ice core from Mt. Waddington, British Columbia (Neff et al., 2012). However, they reported no significant alteration of chemistry from the melt above this layer and no additional analysis of this layer was discussed (Neff et al., 2012). In 2015, a PFA was found on Holtedahlfonna icefield in Northwest Svalbard (Christianson et 260 al., 2015), and in 2019 a PFA at Lomonosovfonna ice field on Svalbard was investigated (Hawrylak and Nilsson, 2019).
According to Pohjola et al. (2002), Lomonosovfonna has high amounts of melt, enough to alter the stable isotopes and chemistry of an ice core taken in 1997, yet it was all refrozen within the annual accumulation layer. Now, 20 years later it is known that there is a PFA on the same icefield (Hawrylak and Nilsson, 2019). Apart from these studies, there have been no other published reports of PFAs on mountain glaciers (Christianson et al, 2015). 265 According to Munneke et al. (2014), firn aquifer formation in Greenland is contingent upon a high amount of accumulation, which helps to insulate the underlying firn from the winter cold wave. Mean annual temperatures in Greenland are well below 0°C and firn aquifers require: (i) latent heat release from meltwater refreezing, to warm the snow and firn to 0°C, along with (ii) meltwater penetration to depths of 10 m or more, to evade the winter cold wave (Munneke et al., 2014). Using 270 the partially preserved section of firn and earlier reports from the IRRP (Grew and Mellor, 1966;Wood, 1963), the estimated accumulation rate at our core site on Kaskawulsh is 1.8 m w.e. per year. This is similar to reported accumulation rates where PFAs have been identified in southeastern Greenland (e.g., Miège et al., 2016).
Temperature records from the previously mentioned weather station 12 km southwest of our core site (Figure 1) at 2600 m 275 a.s.l provide data for summer positive degree days (PDD), which can be taken as a proxy for summer melt extent (Table 2), using the standard approach (Cuffey and Paterson, 2010). These data indicate that there was ample melting every summer from 2014-2018, but with considerable interannual variability. Temperatures in summer 2015 were the warmest of the period 2014-2018, and it is possible that there was deep meltwater penetration and aquifer formation or recharge during that summer. However, we do not know how old the deep aquifer is, or if it experiences annual recharge. It is possible that all of 280 the meltwater refreezes near the surface in cool summers (e.g., 2014), and deep percolation only occurs in years with high amounts of melt (e.g., 2015, 2018).

Density of the upper Kaskawulsh Glacier
In the upper accumulation zone of Kaskawulsh Glacier the density of the upper ~17 m of firn has an average value of 608 ± 2 kg/m 3 in our two cores. The firn extends down at least 36 m from the snow surface, and potentially continues to further 285 depths, giving firn depths of at least 32 m. The average firn density over the full 32 m of the core was 670 ± 2 kg/m 3 . In https://doi.org/10.5194/tc-2020-119 Preprint. Discussion started: 25 May 2020 c Author(s) 2020. CC BY 4.0 License. , Foy (2009 measured the average density of the first 8 meters as 510 kg/m 3 with a standard deviation of 110 kg/m 3 , although this estimate includes some snow. Grew and Mellor (1966) discuss a 15-m core (including snow) that was retrieved at IRRP Divide Camp A. The first 15 m are much less dense than the data we collected, with an estimated average density of roughly 500 kg/m 3 . Their firn densities do not reach 600 kg/m 3 until depths of more than 11 m, while we recorded values in 290 excess of 600 kg/m 3 within the first 5 m. This also implies limited ice content in the IRRP core. One notable ice lens and one ice gland were reported in the Grew and Mellor study (1966). This study also suggests that there may be an indication of meltwater percolation into previous years' accumulation. A recent study identifies the St. Elias Mountains as a region that has warmed in recent decades (Williamson et al., 2020). A downscaled North American Regional Reanalysis (NARR) temperature dataset validated against regional meteorological data from 1979-2016, and shows a warming temperature trend 295 of 0.021 ± 0.012 °C a -1 for elevations between 2500-3000 m (Williamson et al., 2020). This amounts to a warming of roughly 1.2°C at our study site since the 1960s, which would increase the quantity of meltwater, refrozen ice, and meltaffected firn, increasing the density of firn in the St. Elias Mountains.
A range of values for firn density has been reported for alpine glaciers, with considerable variability (Table 3)

. Kaskawulsh 300
Glacier firn is similar to the firn reported in Arctic Canada of 560 kg/m 3 (Zdanowicz et al., 2012) and that of Patagonia (Matsuoka and Naruse, 1999;Shiraiwa et al., 2002). It is more dense than the firn of Svalbard (Pälli et al, 2002;Nuth et al., 2010), but less dense than firn reported from the European Alps and the Himalayas (Table 3). The lower latitude glaciers are likely experiencing a substantial amount of meltwater percolation and refreezing, making ice layers more prevalent than the firn itself. The glaciers in central Svalbard have a different accumulation regime and may have a lower density snow to start 305 with, due to the cold temperatures and low precipitation (depending on the location) (Pälli et al, 2002;Nuth et al., 2010;Eckerstorfer and Christiansen, 2011). Kaskawulsh Glacier firn is similar to the firn at other high latitudes that have not begun to melt significantly, yet.

Implications for geodetic mass balance
The effects of meltwater storage through refreezing or liquid retention on high mountain glaciers will complicate mass 310 balance measurements. Geodetic mass balance measurements are compromised by climate change-induced densification that causes surface lowering of the accumulation zone (Reeh, 2008;Huss, 2013). Mass balance studies in Greenland indicate that changing melt regimes and the subsequent meltwater storage and unknown density of snow and firn pose the greatest uncertainty when modeling the mass balance of large ice sheets (Lenaerts et al., 2019). Research in Greenland proposes that ice-layer formation and the presence of firn aquifers may delay surface run-off until a certain amount of ice layers are 315 present. If ice becomes too extensive, it can form an impermeable barrier that leads to enhanced surface runoff (MacFerrin et al., 2019). These phenomena and effects are not limited to Greenland and the high Arctic. This study demonstrates that Kaskawulsh Glacier also experiences meltwater storage in the form of ice layers and liquid water retention. The firn https://doi.org/10.5194/tc-2020-119 Preprint. Discussion started: 25 May 2020 c Author(s) 2020. CC BY 4.0 License. stratigraphy found on Kaskawulsh Glacier is similar to high accumulation Arctic regions (Pohjola et al., 2002;De La Peña et al., 2012;Bezeau et al., 2013). 320 The accumulation zone of Kaskawulsh Glacier is estimated to have experienced 1.3 ± 0.8 m of surface lowering due to internal refreezing over the period represented by Core 1, which we estimate to be 13 years. This is a conservative estimate, because neither the meltwater retention due to the infiltration ice nor the presence of the firn aquifer is included in this estimate. Firn that is affected by these processes will be denser, with greater surface lowering associated with firn 325 densification, and no change in mass only change in volume. Assuming 13 years of net accumulation (between 2005 and 2018) represented by Core 1, the average thinning rate due to ice-layer formation is 0.10 m/yr, approximately 2% of annual net accumulation. We suspect that the result is representative of the upper accumulation area of Kaskawulsh Glacier and divide region, an area of ~81 ± 12 km 2 or 10% of the accumulation zone area (811 km 2 ; Foy et al., 2011). Additional data is needed to assess the regional applicability, however given that temperature is the primary control on meltwater extent and 330 refreezing (Samimi and Marshall, 2017;Parry et al., 2007), this may be broadly representative of the upper accumulation and Divide region.
Previous research has focused on elevation changes on Kaskawulsh Glacier in recent decades. Foy et al., (2011) found that

340
Over long periods this may be insignificant due to the natural variability in the system and the uncertainties associated with the estimates. For example, Bamber and Riveria (2007) suggest that if the surface lowering due to densification is on the order of centimeters per year then it will not significantly affect geodetic mass balance techniques, but if it is larger then it needs to be taken into consideration if the time period is less than a decade. Additionally, data on short time periods (tens of years or less) may be skewed if there are a few years of significant surface lowering (Huss, 2013;Schneider and Jansson, 345 2004), or if there is a large seasonal variation in glacier surface height, such as southeast Alaska (Pelto et al., 2013). Therefore, the magnitude of surface lowering and the time period of investigation need to be taken into consideration together. Huss (2013) concludes that the uncertainties in densities when converting to volume change can account for 2-15% of mass change. This indicates that there is a substantial need for in-situ density measurements of glacier firn around the world, so that glacier mass balance can be accurately assessed. Estimates of the firn density of Kaskawulsh Glacier are 350 important for geodetic mass balance measurements in the St. Elias region.
This study advances the understanding of the density of firn in this region to better inform geodetic mass balance studies.
The presence of the PFA coincides with findings in high Arctic regions with warming firn. As the climate changes, and with 360 expected increases in temperature and changes in precipitation patterns, the effects of meltwater percolation and refreezing are expected to become more widespread in high-latitude and high-altitude glaciers. The St. Elias Mountains are home to several large valley glaciers (>40 km in length) and several icefields (e.g., Bagley, St. Elias). Mass balance studies are needed for these glaciers and icefields to assess their contribution to sea level rise and changing downstream hydrological conditions, but their size necessitates geodetic approaches for this research. An estimation of the density of firn in this area 365 will aid in geodetic mass balance studies. With one side of the mountain range experiencing a maritime climate and the other side a continental climate, it is essential to consider firn densification processes and the effect of meltwater percolation and refreezing. The firn characteristics at Kaskawulsh Glacier are likely representative of many melt-affected glaciers in such environments.

370
The Kaskawulsh Glacier PFA needs to be more widely studied. The spatial extent and depth of the aquifer is not yet known.
Ground penetrating radar may be a method to investigate the spatial extent of the feature. Use of an electrothermal drill that can drill through water-saturated firn may allow estimations of the depth of the firn aquifer as well as subsequent studies on the potential flow of the water within the aquifer. This region will likely continue to experience increasing amounts of surface melt and refreezing within the snowpack and firn, therefore there is urgency to obtain climate records from this 375 region.
Data availability. Datasets are being submitted to the Scholars Portal Dataverse. Inquiries can be directed to the corresponding author.

380
Author contribution. NO and AC collected field data. AC ran ion analyses, supervised the field campaign and helped with figures. SM contributed to the design, conception, and funding of the study. SM and BM provided supervision during the project. LC provided weather station data. NO analysed the data. NO wrote the manuscript, to which all co-authors contributed.