Grain-size evolution of polar firn: a new empirical grain growth parameterization based on X-ray microcomputer tomography measurements

Abstract

Firn microstructure properties from six different sites in Greenland and Antarctica are investigated by means of X-ray microcomputer tomography. The optical effective radius is calculated from the specific surface area (SSA) and used as a measure of grain size. It is shown that the recently introduced spherical approximation of firn grains using the effective radius R_eff is representative of grain size in the microwave frequency region. The measured profiles show the well-known increase of grain size with depth at all sites, where the increase is largest at near-surface depths. A large variability in grain size on the decimeter-to-centimeter scale as a result of different grain properties of single layers is superimposed on the overall trend at each site. A simple empirical parameterization of grain-size evolution is developed which allows the rapid grain growth in the uppermost layers of the firn to be predicted. The growth is driven by strong seasonal and diurnal temperature gradients. The model can be used to simulate grain-size profiles required by models of firn/microwave interaction (e.g. for retrieval of accumulation rates from satellite microwave sensors) in a more realistic fashion.

1. Introduction

Microwave remote sensing has become an important tool to observe changes in mass balance of the polar ice sheets (Flach and others, 2005; Lacroix and others, 2008). While representative point measurements on the ground over large areas are often lacking, remote-sensing methods can in principle be used for spatial upscaling (Domine and others, 2008). However, for analysis of the resulting images it is necessary to consider microwave/firn interactions in the upper layers of the snowpack, as penetration depths of microwaves into the firn range from cm to ˜100 m (Legresy and Remy, 1998). For this application, knowledge of the grain size plays a crucial role, since the backscattering behaviour (for active microwave imaging) as well as the emissivity (for passive microwave imaging) is strongly influenced by grain size (e.g. Rott and others, 1993; West, 1994). For optical remote-sensing techniques, the surface grain size is equally important, since it influences the albedo (Aoki and others, 2000). However, using microwave sensors, the grain-size evolution with depth needs to be considered. Thus, in order to interpret microwave images, the grain- size profile as well as the interaction of microwaves with the grains has to be considered (Brucker and others, 2010). Therefore, it is important to validate grain-size parameters used in remote-sensing applications with ground truth.

Microstructure measurements of the upper firn column of the polar ice sheets are difficult and time-consuming to obtain. Grain-size measurements were conducted at single locations in Antarctica and Greenland by, for example, Gow (1969), Alley and others (1982), Gay and others (2002), Freitag and others (2004), Gow and others (2004) and Courville and others (2007). Few studies exist that compare firn microstructure at different sites (e.g. Taylor, 1971; Nishimura and Maeno, 1985; Shiraiwa and others, 1996; Gay and others, 2002). The difficulties in comparing the published data lie in the application of different methods of grain-size definition and determination, and the different spatial resolution of the datasets. Another problem is that measurements were taken at different times of the year. Since snow grains show a large variety in their shapes, a number of different parameters are in use to quantify the grain size. Recently, the effective radius determined from the specific surface area (SSA) was introduced as a suitable parameter for grain-size quantification, as it takes the interaction at the air/ice interface in chemical or physical processes into account (Domine and others, 2008). Several alpine snow measurements and experiments were conducted using the effective radius, the SSA and other physical properties as comparable parameters for different snow types and their change with time (Flin and others, 2004; Schneebeli and Sokratov, 2004; Kaempfer and Schneebeli, 2007). SSA profiles of the upper meter of the firn column were measured in Antarctica by Gallet and others (2009, 2011). To date, the most comprehensive overview of available snow properties is given by Domine and others (2008) for alpine snow. An extensive study conducted by Gay and others (2002) for Antarctic snow uses a different grain-size definition, and sampling depths are restricted to the upper 2 m of the snowpack. A comparable and systematic study of polar snow and firn, investigating the grain size in terms of evolution with depth and time under different local climate conditions, does not yet exist. Computer tomography has established itself as a valuable tool for examining snow and firn microstructure (Schneebeli and Sokratov, 2004; Kaempfer and Schneebeli, 2007), being well suited to realistically reproducing snow parameters such as density (Freitag and others, 2004). Using this method, we analyze a set of grain-size measurements from six climatologically heterogeneous sites on the polar ice sheets. We use the effective radius obtained from the SSA to investigate the grain-size evolution as a function of depth and time. From the new dataset, a simple grain growth model is developed which is able to capture the rapid grain growth in the uppermost firn layers for a wide range of polar climate conditions. This model can be used to verify and improve grain-size modeling for remote-sensing applications.

2. Methods

2.1. Firn-core locations

Six surface firn cores, one from Greenland and five from Antarctica, are analyzed. They cover a broad range of annual mean surface temperatures, accumulation rates, elevations and distances to the coast (Fig. 1; Table 1). Annual mean temperatures and accumulation rates are obtained from published data, except for the Antarctic site Depot 700 (DP7; Table 1). The B26 firn core shows intermediate annual surface temperature conditions on the Greenland plateau and medium accumulation rates. In our dataset from Antarctica the location of the FireTrack (FT) firn core from Dome C has the lowest temperatures and accumulation rates. Unfortunately, the available data start from 6.6 m below the surface to larger depths. The DP7 firn core shows a similar annual mean temperature (estimated from Moderate Resolution Imaging Spectroradiometer (MODIS) surface temperature data (LP DAAC, 2009)), and a slightly higher accumulation rate estimated from interpolation of accumulation rates of nearby sites (Isaksson and others, 1999). From Kohnen station (EDML drilling site) two firn cores, B35 and B36, are available, representative of medium-range accumulation rates and temperatures on the Antarctic plateau. While the Hercules Dome firn core represents medium to low temperatures and medium accumulation rates, the B38 firn core from the coastal area shows a remarkably high accumulation rate and a higher mean temperature. With a temperature range between –18° C and –54.5°C and accumulation rates between 0.025 and 1.25 mw.e., our dataset covers a large fraction of climate conditions of polar ice sheets.

Fig. 1. The location of the firn cores.

Table 1. Firn-core locations and local climate conditions

2.2. Density measurements

Density was measured with a vertical resolution of 1 mm using a non-destructive logging system including a Loffel densimeter (Wilhelms, 1996). The radiation source was¹³⁷Cs. Using the measured gamma-ray signal intensity, I, in relation to the signal’s intensity in air, I (0), the mass absorption coefficient µ_ice = 0.085645±0.01 m² kg-¹ (Wilhelms, 1996, 2000) and the diameter d of the medium, the density can be calculated. Details are given by Wilhelms (1996, 2000) and Horhold and others (2011). All measurements were conducted at temperatures of —20°C in the cold laboratory of the Alfred Wegener Institute (AWI), Bremerhaven, Germany. High-resolution density measurements were conducted for all but the Hercules Dome firn core. For each core, an approximate age profile was calculated using the water equivalent depth derived from the density measurements and the known accumulation rate.

2.3 Computer tomography and image analysis

Microstructure images were generated by microcomputer tomography (CT) using a SkyScan CT1072 operating at 80 kV for B38 and a CT1074 operating at 40 kV for all the other firn cores. Cylindric snow samples of 2.5 cm length and 2 cm in diameter (1 cm for B38) were placed on a table in front of an X-ray source. The table was rotated in steps of 0.9°, and each time a shadow image was taken. The resulting set of images was converted to a stack of horizontal gray- value images at 15.73 |am resolution (B38) and 40 |am (all other cores) using a back-projection procedure. For cores B35, B36, B38 and DP7 the uppermost 1–2 m were sampled continuously Below, ˜16 samples were taken on a length of 40 cm for each 1 m interval (Fig. 2). B26 is sampled every meter from the surface. The Hercules Dome dataset contains discontinuous samples with non-equal steps throughout the firn core. The lowermost meter of both B38 and B26 was sampled continuously again. The number of samples in each interval varies from 16 to > 20 (Table 2).

Fig. 2. Firn-core sampling schematic

Table 2. Depth interval sampled, approximate age interval covered by the samples and number of samples at each core site

For image analysis the softwareMAVI (Modular Algorithms for Volume Images) developed by the Fraunhofer Institute for Mathematics, was used (Armbrecht and Sych, 2005). The input data consist of an image stack comprising at least 400 images. From each measured sample a cube of 400 × 400 × 400 voxels (6.3× 6.3 × 6.3mm3 for B38, and 16× 16× 16mm3 for the other cores) was extracted. This size is sufficiently large to be a representative volume for the firn properties considered (Coléou and others, 2001; Kaempfer and Schneebeli, 2007). The datasets were smoothed using a 3 × 3 × 3 median filter. Next, two Gaussian distribution functions were fitted to the gray-value distribution of three images of each stack to determine the threshold value for the segmentation of pore and ice phase. The arithmetic mean of the maximum value of each of the Gaussian distributions determined over three images was used as the threshold gray value for segmentation. Subsequently, an object filter was applied to remove noise, i.e. voxels adding < 1% to the ice or pore phase. All microstructure parameters obtained by MAVI represent the volume structure of each firn sample rather than, for example, single grains.

Porosity was calculated as the relation between pore fraction and total volume of the sample. Density can be computed from the ice fraction times the density of ice (ρ _ice =0.917 g cm ^{− 3}). The chord length l is defined as the average particle or pore extent in the direction of the sample cube’s axes. The SSA represents the ice/air interface per unit mass:

(1)

The SSA is given in cm² g ^{− 1}, S _d denotes the surface density, i.e. the ratio of the total surface to the total volume of the sample, and ρ is the density of the sample (Domine and others, 2008). The so-called effective radius R _eff can be defined as the radius of equivalent-sized spheres with the same SSA. The relation between SSA and R _eff is

(2)

Here it is assumed that the ice phase of the snow and firn can be represented by spheres with a radius r sphere, such that r _sphere = R _eff. To validate this assumption, we used the measured chord length l to calculate the equivalent radius of a sphere:

(3)

The SSA_sphere is calculated by substituting the R eff in Eqn (2) with Eqn (3):

(4)

It reproduces the measured values well, especially given that SSA and chord length are obtained independently (Fig. 3). This supports the assumption that firn microstructure can be realistically represented by spheres with an effective radius R _eff obtained from SSA. Accordingly, in the following we denote grain size by R _eff.

Fig. 3. The SSA obtained from measured surface and density, with the measured chord length of the different sites. The black line is the SSAsphere, calculated from the measured chord length under the assumption of equivalent spheres (Eqn (4)). Note the logarithmic scale of the SSA axis.

3. Measurement Results

3.1. Density

The density of all firn cores except HD was measured using the method described in Section 2.2. For the HD core the density was deduced from the CT-image analysis as described in Section 2.3. The density profiles (smoothed by a running mean of 20 mm) of the firn cores show large fluctuations due to the firn layering (Fig. 4) and differ in the increase of density with depth. B38 shows the highest density, followed by Hercules Dome. B26 starts with lower density, but its density gradient is larger than those of B36 and DP7. FT shows the lowest density.

Fig. 4. Density profiles of the measured firn cores.

3.2. SSA and grain size

The smallest grain sizes are observed at B38, while DP7 and B26 show the largest grain sizes. At greater depths, FT has the largest grain sizes (Fig. 5a). The SSA decreases rapidly in the uppermost meters. B38 shows the highest values at the surface, and DP7 and B26 the lowest (Fig. 5b) relative to each other. At greater depths, FT shows the lowest SSA values. For all sites we find a rapid decrease in SSA in the upper few meters (Fig. 5b). From visual inspection it seems that the strongest increase in grain size with depth in the uppermost meter occurs at B26 and DP7.

Fig. 5. Variation of the SSA and of R eff with depth at the different sites.

4. The Grain Growth Model

4.1. Parameterization of grain growth

Various methods for modeling grain growth processes have been described in the literature. Colbeck (1983), Jordan (1991) and Baunach and others (2001) modeled temperature gradient growth as a process driven by the saturation vapor-pressure gradient. In this case, grains grow only when there is a temperature gradient present within the snowpack, neglecting equi-temperature growth processes. Flanner and Zender (2006), among others, have created complex three-dimensional models of firn metamorphism that consider physical processes of grain growth in great detail but require substantial computation time. For a number of applications, though, it is desirable to use a simple, fast and empirical one-dimensional approach to model firn properties. Maeno and Ebinuma (1983) found that snow grain growth can be described as a pressure sintering process, and there are many examples in the literature where the Arrhenius equation is used to model equi-temperature grain growth as a function of depth (e.g. Alley and others, 1982; Gow and others, 2004; Flach and others, 2005). However, due to a lack of comparable grain-size data, it has been difficult to assess the accuracy of grain-size profiles simulated using the methods mentioned above. In our approach, grain-size evolution is parameterized using a simple steady-state particle growth model. Snow accumulation has been found to vary on the order of tens of km (Eisen, 2008). The sensor footprints of passive microwave radiometers vary between 10 and a few tens of km. In such cases, short-distance variations of the firn layer properties cannot be resolved by remote-sensing measurements. In addition, very few field data are available that are acquired on a spatially dense grid to accountfor local lateral variations of firn layer properties. For remote-sensing applications we assume a homogeneous stratigraphy on horizontal scales of the sensor footprint size and approximate vertical variations of firn properties by mean profiles.

We initialize our firn profile by distinct horizontal layers corresponding to half-year intervals. Firn properties change from layer to layer but remain constant within each layer. The radius r (t) of a particle is then determined from an initial radius r₀ and the growth rate K as a function of time t

(5)

The Arrhenius equation is used to calculate the growth rate K, which is a function of rate constant K 0, activation energy E, gas constant R and absolute temperature T:

(6)

The above parameterization assumes linear particle growth under isothermal conditions. In the literature, values of E vary between 47.0x 10³Jmol^-1 (Gow, 1969) and 42.4x10³ Jmol^-1 (Paterson, 1994). Flach and others (2005) give a K₀ value of 6.75x 10⁷ mm² a^-1, while the K₀ given by Arthern and others (2010) is on the order of 4.1 x10⁶ mm² a^-1. Growth rates calculated from the above range of parameters vary considerably. Budd and Jacka (1989) assumed E and K₀ are temperature-dependent, and Jacka and Li (1994) supplied parameters to empirically fit activation energy E and growth rate K to temperature. The initial value r₀ is usually fixed at an arbitrary value, due to a lack of reference values obtained from the field. We have chosen to empirically fit K to temperature and to develop a parameterization for the initial grain size.

4.2. Firn temperature

From the measured grain-size profiles we aim to find a parameterization of grain growth as a function of mean annual surface temperature T, maximum temperature amplitude at the surface ΔT and mean accumulation rate A. To determine mean grain size dependent on depth in the firn, the measured effective radii are averaged over the measurement intervals of ˜40cm (see Section 2.3; Fig. 2). Firn age t (z) and depth z are related by the velocity v(z) = Aρ_ice/ρ(z) at which a snow layer is buried below the surface (Munk and others, 2003):

(7)

Temperature propagation into the snowpack T (z) is modeled as an exponentially decaying oscillation as described by Paterson (1994), depending on annual mean temperature , amplitude of the seasonal temperature signal at the surface ΔT, thermal diffusivity of snow k, and frequency ω = 2n and phase φ = 0.279 of the seasonal signal. The growth rate depends on the temperature gradient and the absolute temperature. Based on our observations, we assume that the growth rate at the maximum absolute temperature dominates the particle growth processes over the seasonal cycle. For this reason we fixed the phase value at the day when the maximum of the seasonal temperature is reached. Since the focus of microwave remote-sensing applications is on timescales of decades, and the measured signal is also influenced by deeper firn layers not affected by recent surface processes, we neglect near-surface temperature variations on shorter than seasonal timescales.

(8)

Annual mean temperature and amplitude ΔT can be obtained from the MODIS Land Surface Temperature data product (LP DAAC, 2009). Thermal diffusivity can be calculated from thermal conductivity. The latter depends on temperature and crystal structure and is not easy to parameterize (Domine and others, 2008). Therefore, density is used as a proxy for grain characteristics. Here we apply the empirical approach developed by Sturm and others (1997).

4.3. Growth rate

Input data for determining the growth rate as a function of measured radius are mean grain-size profiles with corres-ponding values for temperature T, maximum temperature amplitude ΔT (derived from the temperature profile in Eqn (8)), measured density and firn age. The growth rate K, which is defined as the change of grain size with time, is evaluated numerically for every depth interval. We compared growth rates determined from our dataset with values published by Paterson (1994), and found our values in good accordance (Fig. 6), albeit systematically larger. A possible explanation for this effect lies in the fact that the Paterson values were determined at greater depths, closer to the firn/ice transition. Additionally, different methods were used to determine grain size, and as a result the above comparison can only be approximate. The temperature- dependent activation energy E and rate constant K 0 are difficult to fit numerically to parameters obtained rom measurements, since they are very sensitive to noise in the data. For this reason, we chose to fit the growth rate K directly to T, and determine coefficients for Eqn (9). From the available data, no correlation between K and A can be inferred, therefore K is modeled to depend solely on temperature:

(9)

Fig. 6. Growth rates derived from profile data and from Paterson (1994).

We obtained the following parameter set to estimate the temperature-dependent growth rate:

(10)

4.4 Surface grain size

For remote-sensing applications, it is convenient to define the initial grain size r ₀ as an average radius over a subsurface layer of predefined thickness, which differs from the grain size at the time of snow deposition and is influenced by the metamorphism processes that took place between the time when the lowest part of the subsurface layer was deposited and the time of the remote-sensing measurements. To quantify the temperature and accumulation dependency of the surface grain size, initial radii r₀ and their standard deviation σ_r0 were estimated from the measured grain sizes of the upper 30 cm interval of each measured individual profile. This was done for all available grain-size profiles except the FT profile since that starts at a depth of 6 m.

A multiple linear regression applied to the Table 3 data yields the following relationship between initial grain radius /o, mean annual temperature (°C) and accumulation rate A (m w.e. a^-1).

(11)

with

(12)

Table 3. Initial radii, radius standard deviation, temperatures, accumulation rates and number of samples N

In order to test the accuracy of our fit, we compared initial radii determined from Eqn (11) with measured values. It can be seen that our multiple regression approach reproduces the measured values with sufficient accuracy (Fig. 7). The grain size of the near-surface layer r₀ seems to depend on temperature and accumulation rate. Larger water-vapor transport and thus faster grain growth are linked to higher temperatures (Domine and others, 2008). Higher accumulation reduces the time the snow grains are subject to a strong temperature gradient and slows the grain growth process.

Fig. 7. Correlation between simulated and measured surface grain sizes.

4.5. Modeled grain-size profiles

Figure 8 shows that the rapid grain growth in the upper layers influenced by a strong temperature gradient is represented more realistically by the new approach. Our model overestimates grain growth for B38, a site with an extremely high accumulation rate of 1.25 mw.e. a^-1 and a comparably high mean annual surface temperature of — 18°C. For DP7, the extreme particle growth is caused by a very low accumulation rate of ˜0.045 mw.e. a^-1 and, in consequence, the long exposure time of snow layers to a large temperature gradient. This effect is underestimated in our model. For intermediate polar climate conditions, the simulated grain-size profiles closely fit the measurements. Figure 8 compares model results with measurements and grain-size parameterizations by Paterson (1994) and Zwally and Li (2002). Since start values for the grain growth were not given in both examples, we calculated initial radii r₀ from Eqn (11) for better comparison. Both approaches also model grain sizes using Eqns (5) and (6), albeit with different values for growth rate and activation energy The parameters given by Paterson (1994) are optimized to reflect particle growth under isothermal conditions (see Fig. 6; Section 4.1) and hence does not represent conditions found in the uppermost firn layers in the presence of a strong temperature gradient. The work by Zwally and Li (2002) results in a temperature-dependent parameter set for growth rate and activation energy and is based on data from field studies and laboratory experiments Jacka and Li (1994) obtained under isothermal conditions. Hence, results from this parameter set are expected to deviate from observations of particle evolution in near-surface layers.

Fig. 8. Modeled grain-size profiles for six polar firn cores

The competing impact of temperature and accumulation rate leads to overestimation of the B38 samples and underestimation of the DP7 and FT samples. This indicates a lower sensitivity of the model parameters to extreme conditions, but the overall agreement of data and model results is good. Possible sources of error in the development of the parameterization stem from uncertainties in accumulation rates and temperatures, which are difficult to quantify Even though the results of the linear regression for the initial radius r₀ (Fig. 7) seem to be nearly perfect, it is difficult to say whether the parameterization approach remains valid when more datasets are considered. A comparison of modeled grain radii with model results from Flanner and Zender (2006) for temperature gradient growth yields comparable grain sizes.

5. Conclusion

From computer-tomography measurements we showed that firn grains can be approximated as spheres for remote- sensing applications. The use of an effective (or optical) radius can be reproduced by different measurement methods and can be used to parameterize grain size. A simple, empirical parameterization of grain growth processes based on the Arrhenius equation is introduced. It incorporates effects of the strong particle growth in the upper layers of the snowpack caused by thermal gradients as well as equi-temperature growth at depths no longer influenced by the seasonal temperature change. The new model includes a parameter set for estimating the surface grain size as a function of accumulation rate and mean air temperature at the surface. Since the firn cores used to derive our model parameter set do not exceed ˜20 m depth, we cannot predict the accuracy of our model results for greater depths. Therefore we recommend this model for remote-sensing applications for which firn/microwave interactions have to be considered in the upper layers of the snowpack down to 2030 m. Because microwave radiation is sensitive to particle size, a more accurately simulated firn profile can be used to improve many applications (e.g. the retrieval of snow accumulation rates or an inversion of surface grain sizes from microwave satellite data). Since we used data from six firn cores from Greenland and Antarctica which represent a heterogeneous set of environmental conditions, we can safely assume our model is applicable to many regions of Greenland and Antarctica. In order to further test the validity of our relation for estimating surface grain sizes, more data will be needed.