Multi-tracer study of gas trapping in an East Antarctic ice core

. We study a ﬁrn and ice core drilled at the new "Lock-In" site in East Antarctica, located 136km away from Concordia station towards Durmont d’Urville. High resolution chemical and physical measurements were performed on the core, with a particular focus on the trapping zone of the ﬁrn where air bubbles are formed. We measured the air content in the ice, closed and open porous volumes in the ﬁrn, ﬁrn density, ﬁrn liquid conductivity and major ion concentrations, as well as methane concentrations in the ice. The closed and open porosity volumes of ﬁrn samples were obtained by the two independent methods 5 of pycnometry and tomography, that yield similar results. The measured increase of the closed porosity with density is used to estimate the air content trapped in the ice with the aid of a simple gas trapping model. Results show a discrepancy, with the model trapping too much air. Experimental errors have been considered but do not explain the discrepancy between the model and the observations. The model and data can be reconciled with the introduction of a reduced compression of the closed porosity compared to the open porosity. Yet, it is not clear if this limited compression of closed pores is the actual mechanism 10 responsible for the low amount of air in the ice. High resolution density measurements reveal the presence of a strong layering, manifesting itself as centimeter scale variations. Despite this heterogeneous stratiﬁcation, all layers, including the ones that are especially dense or less dense compared to their surroundings, display similar pore morphology and closed porosity as function of density. This implies that all layers close in a similar way, even though some close in advance or later compared to the bulk ﬁrn. Investigation of the chemistry data suggests that in the trapping zone, the observed stratiﬁcation is partly related to the 15 presence of chemical impurities. system was restored by OM and LA, and the pycnometry measurements were performed by XF and KF. The tomography scans and segmentations were carried out by HL and MS. The air content measurements were performed by 20 VL. The Lock-In project was supervised by PM and XF. The code development for data processing and modeling was performed by KF and PM. All authors contributed to the discussion and interpretation of the data. The article was written by KF with the help of all co-authors.

Among the polar regions, East Antarctica is of particular interest for ice core drilling as this is where the oldest ice is retrieved. While several studies focus on the cold and arid sites of this region, none of them to our knowledge attempts to encompass a detailed description of stratification, closed porosity and gas trapping within the same firn. As particular 5 examples, the study of Fujita et al. (2016) describes the link between density variability and ionic content in three firns core near Dome Fuji without addressing pore closure, while the studies of Schaller et al. (2017) and Burr et al. (2018) propose detailed descriptions of the pore network closure but without discussing stratification.
Our study thus aims to provide a detailed description of chemical and physical properties along the trapping zone of an Antarctic firn core. This core was drilled at the site of "Lock-In" in East Antarctica. It was chosen as its temperature is similar 10 to the cold sites of Dome C and Vostok where deep ice cores have been drilled. Nonetheless, "Lock-In" displays a higher annual precipitation rate than Dome C and Vostok. When possible, the measurements were performed with a resolution allowing to observe stratification and variabilities at the centimeter scale in the firn. Using high-resolution datasets, we aim at investigating the impact of firn stratification on gas trapping.

Lock-In study site
The ice core studied in this article was drilled in January 2016 at a site called "Lock-In", located on the East Antarctic plateau, 136 km away from Dome C along the traverse road joining Concordia and Dumont d'Urville stations (coordinates 74 • 08.310 S, 126 • 09.510 E). As an East Antarctic plateau site, "Lock-In" exhibits a fairly low mean annual temperature of −53.15 • C (measured at 20 m depth in the borehole), slightly higher than at the deep drilling sites of Dome C and Vostok. The site has 20 an elevation of 3209 m above sea level, and the surface atmospheric pressure was estimated to be 645 mbar using the value at Dome C of 643 mbar (data from the Automatic Weather Station Project) and a pressure elevation gradient of 0.084 mbar m −1 (evaluated using pressure differences between D80 and Dome C). "Lock-In" is also characterized by a relatively low accumulation rate of 3.6 cm.we.yr −1 . The accumulation rate was evaluated by locating volcanic events with solid conductivity measurements. Hence, "Lock-In" has a similar temperature but a higher accumulation rate than Vostok and Dome C. 25 The drilling operation retrieved ice down to 200 m below the surface. The amount of material was sufficient in the 10 cm diameter core to perform several parallel measurements, which are described below. During the drilling procedure, air was sampled from the firn open porosity along the firn column. The last pumping was performed at 108 m depth, giving an indication of where the firn reaches full closure and transforms into airtight ice. Visual observation of the firn core did not reveal any obvious ice layer.

High resolution density profile
Density is a parameter of particular interest when studying the enclosure of gas in polar ice. Indeed, the densification speed 5 measures the amount of compaction of the firn and the reduction in pore volume. Hence, the density was continuously measured in the "Lock-In" firn and ice core. The measurements were performed on the whole 10 cm diameter core from 6 m down to 131 m below the surface, at the Alfred Wegener Institute, Bremerhaven, Germany. The measurement technique is based on the absorption of gamma rays by the ice phase while traversing the core. The density is then derived from the ray attenuation using Beer's law. The method of measurement is described in details in Hörhold et al. (2011). 10 To be applicable, Beer's law requires the estimation of the core diameter d and the ice specific absorption coefficient µ ice .
The diameter of the core was measured using a caliper for one meter long sections. However, due irregularities this value might not be representative of the entire section. Therefore, ice core parts where ice was missing (for instance parts broken during drilling) were removed from the final dataset. Finally, the value of the specific absorption coefficient µ ice was calibrated for each new scan using cylinders of pure ice with known diameters. One should note that the results are quite sensitive to 15 the determination of d and µ ice . Thus, the high-resolution density data will not be discussed in terms of absolute values but rather in terms of centimeter scale variations, which are better constrained. This measurement method allows to retrieve density variations at the sub-centimeter scale.

Pycnometry measurements
Pycnometry is an experimental method allowing the measurement of closed and open porosity volumes in a firn sample (Stauf-20 fer et al., 1985;Schwander et al., 1993). In this context, a pore is determined as open if it reaches the edge of the sample. The pycnometry experimental set-up is composed of two sealed chambers of known volumes V 1 and V 2 , with a pressure gauge connected to the first chamber and a valve allowing to isolate V 2 from V 1 . The sample is placed in chamber V 1 , while chamber V 2 is isolated and vacuum pumped. Then, the two chambers are connected, resulting in a homogeneous pressure in both chambers lower than the initial pressure in V 1 . This pressure drop is linked to V s , the volume rendered inaccessible to gases by the where R = P /P with P and P being the pressure value in chamber V 1 respectively after and before the dilation. The application of the equation requires the knowledge of the volumes V 1 and V 2 . The two volumes V 1 and V 2 were calibrated using stainless steel balls and bubble free ice pieces of known volumes. The derivation of Equation 1 and the description of the 5 calibration procedure are provided the Supplement Sections S1.1 and S1.2.
Then the closed and open porosity volumes are computed given that: where V c is the closed pore volume, V o is the open pore volume, V ice is the volume occupied by the ice phase, V cyl is the 10 volume of the sample, M is the mass of the sample, and ρ ice is the density of pure ice. The density of pure ice was taken as 0.918 g cm −3 , for a recorded temperature of −8 • C during the measurements (Bader, 1964;Goujon et al., 2003).
A total of 80 samples were taken within the "Lock-In" core, with depths ranging from 76 m to 130 m. The measured samples were cylinders with typical diameters and heights of 4 cm, obtained with a vertical lathe at the Institut des Géosciences de 15 l'Environnement (IGE), Grenoble, France. The radius and height were measured for each sample using a digital caliper in order to estimate the volume of the firn sample V cyl . To estimate the ice volumes V ice , the samples were weighted. Finally, the densities ρ of the samples were computed as ρ = M/V cyl .
The treatment of uncertainties is reported in the supplementary material Section S1.3. One of the main results of this uncertainty analysis is the large uncertainties associated with the sample volumes V cyl due to their sensitivity to the radii measurements. 20 The large uncertainty on V cyl affects the determination of the density ρ with a typical uncertainty of 0.01 g.cm −3 .
Following the definition of an open pore from a pycnometry point of view, some pores that are closed in the firn (they no longer connect to the atmosphere) might appear open in small samples. For instance, the blue pore in Figure 1 is closed when taken in the whole firn column, but open when observed in a small cylindrical sample. This opening of closed pores when sampling Figure 1. Depiction of a cylindrical firn sample, numerical slices in tomographic images, and pores. The yellow and orange slices, represent two overlapping slices. The blue pore illustrates the cut-bubble effect: a pore initially closed in the firn column appears open in the cylindrical sample. The green pore is fully enclosed in the cylindrical sample, and is thus considered closed in the slices even if it intersects with their boundaries.
firn is known as the cut-bubble effect, and will be discussed in section 3.1.2.

Tomography measurements
Computed X-ray tomography is a technique gaining traction for the study of firn densification and pore closure (Barnola et al., 2004;Freitag et al., 2004;Gregory et al., 2014;Schaller et al., 2017;Burr et al., 2018). Briefly, several 2D X-ray scans of the 5 samples are performed under different angles. Then, based on the contrast of X-ray absorption between the ice and pore phases, tomography techniques are able to reconstruct a voxelized 3D model of the scanned samples. Ten samples from the "Lock- In" firn core were measured at the Institut für Schnee und Lawinenforschung (SLF), Davos, Switzerland, using a commercial tomography apparatus. The tomography samples have the same geometry as the pycnometry ones, and were also prepared at IGE using a lathe. Four of the samples measured using tomography were also measured by pycnometry. The obtained 10 3D models are gray scale images, composed of square voxels of 25 µm in length. The gray scale images were automatically segmented by fitting the gray value histogram to three Gaussian curves (based on Hagenmuller et al., 2013), producing binary images of pore and ice voxels. Using the 3D images, one can then compute the relative density to pure ice ρ R = ρ/ρ ice of the samples as the ratio of ice voxels to the total number of voxels. The relative density represents the volume fraction occupied by the ice phase in the firn sample. Throughout the article, we use relative densities instead of densities (mass to volume ratio).
This was chosen as relative densities are not sensitive to temperature, and thus allow an easier comparison between sites.
Tomographic images give access to the pore network, and can thus be used to determine the volume of closed and open porosity.
As with the pycnometry measurements, a pore is considered open if it reaches the boundary of the sample. Using the software ImageJ and the plug-in Analysis_3D (Boulos et al., 2012), we differentiated the open and closed pores in the 3D images. The volume of the open and closed porosity is then simply computed by summing the volume of the voxels belonging to an open 5 or closed pore.

Air content
Air content, denoted AC in this article, is a measure of the quantity of air enclosed in polar ice. It is expressed in cm 3 of trapped air at standard temperature and pressure (STP) per gram of ice. Air content is directly related to the volume of pores during closure, and can thus be used to estimate the density at which a particular firn layer closed. Air content was measured in the 10 "Lock-In" firn core using the method described in Lipenkov et al. (1995). In order to correct the AC values for the cut-bubble effect, an estimation of the re-opened volume has been performed (based on bubble size measurements, Martinerie et al., 1990).
In total 10 samples were measured, including 4 replicates. Six of these samples were taken around 122 m depth, with a specific focus on stratification, and the four others were taken deeper at 145 m.

15
Several studies highlight the link between density heterogeneities and ionic content in firn layers (Hörhold et al., 2012;Freitag et al., 2013;Fujita et al., 2016). To study the presence of this link in the trapping zone of the "Lock-In" firn, chemical analyses were performed on this part of the core. The measurements were performed at the British Antarctic Survey, Cambridge, United-Kingdom, using a continuous flow analysis (CFA) system. Sticks were cut from the center of the core (34 × 34 mm), covering depths from 90 to 115 m. The sticks were then continuously melted, and the meltwater analyzed for conductivity. Additional 20 meltwater was collected in vials, using a fraction collector, for later analysis with ion chromatography (IC). Four one meter sections were analyzed at 2.5 cm resolution using a dionex reagent-free ion chromatography system in a class-100 cleanroom.
The measured species, associated resolutions, and depth coverages are summarized in Table 1.

Continuous methane measurements
During the drilling operation, about 100 m of airtight ice was retrieved. This ice was analyzed for methane concentration in 25 enclosed bubbles, using the continuous flow analysis system of IGE coupled with a laser spectrometer SARA based on optical- feedback cavity enhanced absorption spectroscopy (OF-CEAS; Morville et al., 2005). Ice cores were cut in 34 by 34 mm sticks, and melted at a mean rate of 3.8 cm min −1 . The description of this system and data processing is available in Fourteau et al.

Dataset synchronization
This study produced three high resolution and continuous datasets: density, chemistry and methane. Although the chemistry data are composed of various species, they are all based on a unique depth scale and were thus treated as a single dataset for synchronization. In order to produce precise depth scales, notes and measurements were taken when beveled breaks had to be 10 removed from the analyzed ice. Still, uncertainties on the depth scales of the order of a centimeter remain. Since this article will discuss centimeter scale variability, it is necessary to synchronize the datasets to properly observe the presence or absence of covariation. The synchronization was performed manually, by allowing depth shifts of the order of one centimeter between core breaks.
3 Results and discussion 15

Closed porosity ratios
Deep firn densification is characterized by the transformation of the open pore network into bubbles encapsulating the air.
Notably, pore closure is quantified through the evolution of the closed porosity ratio CP , defined as the ratio between the closed pores volume and the total pore volume (Mitchell et al., 2015;Burr et al., 2018). This closed porosity ratio falls to zero when the porosity is fully open and reaches 1 when the firn sample porosity is fully closed, making it a simple indicator of pore 20 closure.
For the pycnometry measurements, CP was simply computed for each cylindrical sample from the measured open and closed volumes. On the other hand, thanks to the explicit representation of the pore network, tomography data can be used to access variations at scales finer than the sample size. Hence, the tomography images were numerically divided into 1 cm thick slices, as depicted in Figure 1. The spacing between two slices was set as 0.5 cm, meaning that slices are overlapping. This overlapping was implemented as an equivalent to a running window analysis. The closed porosity ratio, as well as density, were then 5 determined for each slice. To minimize the impact of the cut-bubble effect, the distinction between closed and open pores was performed on the full cylindrical images as displayed in Figure 1. Moreover, the top and bottom slices of each sample were discarded since, as in described Section S1.4 of the Supplement, they are the most affected by cut bubbles. The values of closed porosity obtained using both pycnometry and tomography are displayed Figure 2-a against relative density. Burr et al. (2018) proposed to use the connectivity index as an indicator of pore closure. This index is defined as the ratio 10 between the volume of the largest pore and the volume of all pores in a given sample, and thus equals one for a fully open sample and sharply drops near zero for partially closed samples. Burr et al. (2018) report that it is less sensitive to the cutbubble effect. However, this index cannot be used in models of gas transport and trapping as it does not explicitly represent the volume of closed pores. Therefore, we will use the closed porosity ratio, rather than the connectivity index as an indicator of closure in this article.

Consistency between pycnometry and tomography
The first closed porosity datasets have been obtained using pycnometry or similar methods (Schwander et al., 1993;Trudinger et al., 1997;Goujon et al., 2003;Mitchell et al., 2015). However, the last decade has seen the development of tomography for the study of firn, and estimations of porosity closure based on tomographic images (Barnola et al., 2004;Gregory et al., 2014;Schaller et al., 2017;Burr et al., 2018). It is thus important to evaluate the consistency between the two types of methods.

20
A comparison of the results obtained for "Lock-In" in Figure 2-a indicates an inconsistency between the two datasets, with the tomography data showing a steeper increase of closed porosity with density. To investigate this difference, we focused our attention on four specific cylindrical samples that have been measured both with the pycnometry and tomography methods.
Comparison of the densities of these samples confirm a disagreement between the two methods, with systematically higher density values obtained with tomography. While the differences fall within the uncertainty range of the pycnometry measure-25 ments, a stochastic error cannot explain the observed systematic discrepancy. To estimate a possible bias of the tomography based density, four gray scale images were manually re-segmented. Relative deviations in density values were of the order of 9.10 −4 , that is to say 0.1%. These deviations are similar to the uncertainties reported by Burr et al. (2018). Moreover, Burr from 30 to 12µm does not change the density by more than 0.1% for relative densities in the range considered here. Therefore, it appears that the densities obtained with tomographic scanning are fairly robust. On the other hand, the densities obtained by size measurements and weighting are sensitive to the volume estimation, particularly through the radius measurement. A close look at the cylindrical samples, notably using the tomography 3D models, revealed rough surfaces and some broken edges.

5
An overestimation of the measured radius by 0.05 mm is enough to explain the density discrepancy between the four samples measured with both pycnometry and tomography. Our understanding is that while measuring with a caliper, the operator tends to overestimate the actual radius of the sample because of the irregular geometry, leading to a density underestimation. We thus corrected the whole pycnometry dataset by considering a systematic overestimation of the radius by 0.05 mm. This is consistent with the observation of sample surfaces using tomography models, reduces the overall discrepancy between the two 10 methods and is within the uncertainty range for radius measurements. This radius correction has an influence on the density values but also on the open pore volumes (Equation 2). The corrected pycnometry results are displayed Figure 2-b, along with the unaltered tomography results. This correction is sufficient to explain the discrepancy in density between the two methods, while also improving the agreement on measured closed porosity ratios. Thus, after correction, pycnometric and tomographic methods produce consistent results of closed porosity.

Cut-bubble effect
As  (Martinerie et al., 1990;Lipenkov et al., 1997;Schaller et al., 2017). The impact of the cut-bubble effect is clearly visible Figure 2-a and b: while closed porosity ratio should reach 100% for high densities (corresponding to a fully closed material), an upper limit of CP = 90% appears in the data. This is due to the presence of cut bubbles at the surface of the cylindrical samples, even in fully closed ice. This cut-bubble effect can be partially mitigated by 10 the use of large samples.
To correct for this effect, one has to estimate the fraction f rac of bubble volume opened during sampling. Then, the true close porosity ratio CP true is given by: perature and 3.0 cm.we.yr −1 accumulation rate). Using tomography, they scanned large firn samples of 10 cm in diameter and 5 cm in height. They then extracted 5 cm diameter cylinders from the large samples. By comparing the state of pores before and after the extraction, they were able to determine which bubbles have been opened. Their findings indicate that up to 60% of the closed porosity gets re-opened near the density of full closure. Moreover, they find that still 30% of bubbles are re-opened when measuring deep ice samples with relative density close to 0.94 (Figure 4 in Schaller et al., 2017). 25 However, when applied to the "Lock-In" closed porosity ratio, the correction proposed by Schaller et al. (2017) leads to unphysical results. For high relative densities, direct application of their correction yields a closed porosity ratio well over 100%.
This indicates an overestimation of the fraction of re-open bubbles for the "Lock-In" firn.
Two methods were used to estimate the fraction of cut bubbles specifically for the "Lock-In" samples. First, for depths well below the firn-ice transition, one knows that the samples are fully closed. Using pycnometry, we thus closed and open pores in the sub-sample was then performed following the method of section 2.4. Then, by comparing with the full sample image, one is able to measure the fraction of closed pores that have been opened. To produce the sub-samples we trimmed the original cylinder to remove 1 cm in height and 0.5 cm in diameter. The goal is to remove most of the immediate 10 boundary effect described Section S1.4 of the Supplement, while conserving a size and a geometry close to the full-sample.
The method would yield more robust results if, similarly to Schaller et al. (2017), larger volumes would be reduced to 4 cm diameter cylinders. However, large volume tomography was not available for the "Lock-In" firn core. The deduced fraction of re-opened bubbles is displayed in Figure 3 as a function of relative density. Our findings suggests that up to 25% of the bubble volume is re-opened in our samples for relative densities closed to 0.90. This value drops to around 7.5% for high densities.

15
In order to correct the CP values obtained by both pycnometry and tomography, we derived a linear piecewise relationship relating the relative density of a sample and the cut bubbles fraction f rac, displayed in black in Figure 3. This linear law was derived to match the estimated re-open fraction, while still resulting in physically sound closed porosity ratios, i.e. a monotonous increase without ratios well above 100 %. We then applied the correction from Equation 3   Unfortunately, there is a lack of available data for relative densities ranging from 0.86 to 0.89, as no tomography sample falls in this range. To test the possible sensitivity to a bad estimation of the re-open fraction at these densities, we derived a second correction law with a re-open fraction set at about 25 % for all the relative densities below 0.91. This correction and the corrected data are displayed in the Supplement Section S1.5. Application of this correction does not significantly change the final closed porosity data, as low closed porosity values are obtained for relative densities below 0.91.

Reconciling air content and closed porosity
The air content measurements in the "Lock-In" core indicate a value of 0.0915 ±0.0017 cm 3 .g −1 of air at standard temperature and pressure per gram of ice based on four samples at 145 m depth, and of 0.0874 ± 0.0013 cm 3 .g −1 based on six samples at core (Lipenkov et al., 1997), but not observed in Dome C ice cores (Raynaud et al., 2007). This suggests a potential change of local conditions of the "Lock-In" and Vostok sites during the late Holocene, that did not occur at Dome C. where V i is the total porous volume at air isolation (expressed in cm 3 per gram of ice), AC is the air content (in cm 3 STP per gram of ice), T i and P i the temperature and pressure of pores at isolation, T 0 = 273.15 K and P 0 = 1013 mbar are the standard temperature and pressure, and ρ c the estimation of the density during closure. Note that V i is not the porous volume when the firn reaches full closure. Indeed, the volumes of the closed pores evolve in the trapping zone due to compression. Therefore, V i is the sum of the individual pores volume when they close, and is larger than the porous volume at complete closure. This 5 means that ρ c is smaller than the density at which full closure is reached, and should rather be seen as the typical density at which bubbles are formed when the air pressure starts to deviate from the atmospheric value. Finally, ρ c can be expressed as a relative density ρ R c by assuming a density of pure ice at −53.15 • C equal to 0.924 g cm −3 (Bader, 1964;Goujon et al., 2003). P0ρ , where R is the ideal gas constant. To apply the model to the "Lock-In" core, the measured closed porosity ratios have been fitted by a spline to obtain a smooth porosity curve (in blue in Figure 4). Moreover, the high resolution density data were fitted by a polynomial function of degree 5 to produce a smooth density curve. Using 20 these smooth closed porosity and density curves in the model, leads to an estimated air content value of 0.102 cm 3 g −1 . This is about 10% higher than the measured values, and well outside the uncertainty range of the measurements.
To explain the discrepancy between data and model, we first investigate the impact of the input uncertainties on the modeled air content value. A first source of discrepancy between the modeled and observed air content values could be an underestima- 25 tion of the densities in the closed porosity-density relationship. Indeed, if densities are underestimated, then the model would trap air at too low densities, resulting in a too large bubble volume and too much air. To test this hypothesis, we shifted the smooth closed porosity ratio curve towards higher densities in order to obtain air content values similar to the ones measured in the ice. To achieve this result, the relative densities have to be shifted by an amount of about 9 . 10 −3 . The resulting closed porosity curve is displayed in orange in Figure 4. However, such a large shift cannot be reconciled with the pycnometry and tomography measurements. Indeed, the absolute errors on relative densities from the tomography apparatus, estimated as 9 . 10 −4 in Section 3.1.1, are much smaller than 9 . 10 −3 .
A second explanation could be a bad estimation of the cut-bubble correction. Indeed, if we overestimated the cut-bubble effect in the middle of the trapping zone, the model would also trap too much air at low densities and thus overestimate the air 5 content. To test this hypothesis, we used the extremal case of a cut-bubble correction corresponding to only 7.5% of opened bubbles throughout the trapping zone. The modeled air content then decreases to 0.096 cm 3 g −1 , still much higher than the measured values of 0.0874 and 0.0915 cm 3 .g −1 . A bad estimation of the cut-bubble correction thus cannot explain the model overestimation.
Finally, the discrepancy could be due to a poor estimation of the surface atmospheric pressure and of the temperature at the 10 "Lock-In" site. Using a surface pressure of 643 mbar, similar to Dome C, and a temperature of −50 • C, the modeled air content only decreases to 0.0994 cm 3 .g −1 . Therefore it appears that the discrepancy between the modeled and measured air content is not due to errors in the input data of the model.  (Stauffer et al., 1985;Martinerie et al., 1992) they argue that the presence of impermeable layers above an open firn layer could retain air in the open porosity. However, this mechanism yields an increase in the amount of air trapped, and can therefore not explain the low value measured in "Lock-In" ice.
Another explanation could be the wrong representation of closed pores compression in the trapping zone. Originally, the model 20 is based on the iso-compression of bubbles with the rest of the firn, that is to say that bubbles compress at the same rate as the open porosity. However during compression, the pressure inside bubbles increases with the diminution in bubble volume. This potentially makes bubbles less easily compressed than the open porosity. To test the sensitivity of the total air content towards bubble compression, we modified the Equation 8 of Rommelaere et al. (1997) to reflect a limited compression of bubbles. For this purpose, we introduce a new parameter α such that when the total porosity decreases of X% during compression, the 25 closed porosity only diminishes by αX%. With this new parameter, the equation governing q b air , the molar quantity of air in bubbles, now writes: where, similarly to Rommelaere et al. (1997), order to obtain the measured air contents of 0.0874 and 0.0915 cm 3 g −1 using the measured porosity data, the α parameter has to be set to 35 and 50 % respectively. This indicates that a smaller rate of compression for the closed pores can reconcile the air content and closed porosity measurements. The optimal α value is sensitive to the cut-bubble correction applied to the closed porosity data. If we use a closed porosity data corrected for a maximal re-opened volume fraction of 20%, the α parameter has to be set to close to 70% to reproduce the measured air content values. Yet, even if a smaller compressibility of bubbles can 20 reproduce the measured air content, it is not entirely satisfactory as it currently lacks supporting observations. The driver of pore compression in the firn is the difference in stresses between the ice and pore phase (Lipenkov et al., 1997). As reported by Martinerie et al. (1992), the difference in driving stresses between open and closed pores compression is lower than 8%, and it is not clear how this difference translates in terms of bubble compressibility. It is therefore possible that another mechanism is responsible for the discrepancy between the modeled and observed air content. The explanation could be related to the release 25 of gases though capillaries not observable with pycnometry or tomography techniques (Huber et al., 2006;Severinghaus and Battle, 2006). Finally, we also wonder if some bubbles might get re-opened in the firn due to the building pressure inside them.
Such a re-opening of bubbles in the firn would release the gas enclosed at low density and lower the air content value. The pycnometry measurements also need to be corrected for the cut-bubble effect. However, contrary to "Lock-In" we do not have tomographic images of the Vostok trapping zone that could be used to quantify the re-opened volumes. Thus, we apply the correction previously derived for "Lock-In" in Section 3.1.2 to the raw Vostok pycnometry data. The application of this 5 correction yields physically sound closed porosity values, with continuously increasing closed porosity ratios up to 1. The uncorrected Vostok data are displayed in the Supplement Figure S7. Finally, the Vostok corrected closed porosity data have been fitted to obtain a smoothed closed porosity curve to be used in the gas trapping model.
Similarly to "Lock-In", a direct application of the original Rommelaere et al. (1997) gas trapping model with iso-compression of bubbles yields a predicted air content of 0.0960 cm 3 g −1 much higher than the measured value of 0.0862 cm 3 g −1 . This 10 provides further insights that the gas trapping process might not be well represented in the Rommelaere et al. (1997) model.
Again, the introduction of a limited compressibility of the bubbles with a coefficient α = 50 % is able to reconcile the modeled and observed air content values.

Comparison with other sites and porosity laws
The comparison in Figure 5 of the measured closed porosity data for "Lock-In" with the proposed parametrization for the B53 site by Schaller et al. (2017) reveals a discrepancy between the two East Antarctic sites. While the B53 closed porosity curve reaches full closure for relative densities of 0.90, the data we obtained for "Lock-In" reach closure for relative densities closer to 0.92. Therefore, we were not able to reproduce one of their main findings, namely that firn should close at a common 5 critical relative density of 0.90, independently of the studied sites. This discrepancy could be due to an underestimation of the cut-bubble effect for the "Lock-In" data set. To test this hypothesis, we derived a new cut-bubble correction that reproduces the closed porosity parametrization proposed by Schaller et al. (2017). As seen in Figure 3, this new correction assumes a very large re-open fraction around relative densities of 0.90. This is consistent with the results of Schaller et al. (2017), but inconsistent with the re-opened fraction of 24% deduced from the tomography image at a relative density of 0.9025. It is therefore 10 possible that the true "Lock-In" closed porosities are similar to the ones observed by Schaller et al. (2017) at B53, but a better estimation of the cut-bubble effect is required to answer this question. However, a closed porosity curve similar to the one observed in B53 would trap a higher quantity of air, and further increases the discrepancy between the modeled and measured air content described in Section 3.2. Even in the extreme case of incompressible bubbles (α = 0), the gas trapping used Section 3.2 predicts an air content value of 0.0974cm 3 g −1 , larger than the values of 0.0874 and 0.0915 cm 3 .g −1 measured in "Lock-In".

15
Based on previous pycnometry measurements, JM Barnola proposed an analytical formulation to represent the closed porosity of firn. This parametrization is fully characterized by the total porous volume V i (or equivalently the density during closure ρ c ) as a single parameter. The derivation of this parametrization is not published in the scientific literature, but Goujon et al. (2003) included this Barnola parametrization in their model (Equation 9 of Goujon et al., 2003) and used the linear relation-20 ship between site temperature and V i by Martinerie et al. (1994) to model closed porosity data using temperature as the only parameter. Application of the Goujon et al. (2003) parametrization using the "Lock-In" site temperature yields a predicted total porous volume V i of 0.110 cm 3 g −1 and the closed porosity curve displayed in orange in Figure 5. The resulting closed porosity curve for "Lock-In" clearly underestimates the closed porosity ratios. In the case of "Lock-In", the two air contents values can be used to estimate the total porous volume V i at 0.109 and 0.114 cm 3 g −1 , indicating a possible underestimation of V i by the 25 temperature relationship. Using the V i value of 0.114 cm 3 g −1 in the Barnola parametrization results in the closed porosities displayed in black in Figure 5. This curve is more consistent with the pycnometry and tomography data, but still displays a clear underestimation of the closed porosity ratio. Two problems thus arise when using the Goujon et al. (2003) parametrization to predict the Lock-In closed porosities. First, the estimation of the total porous volume V i using the temperature relationship might be underestimated, which in turns tends to predict pore closure at too high densities. Then, even with a better constrained V i , the Barnola curve appears to underestimate the closed porosity ratios in the "Lock-In" firn.
Finally, the closed porosity data measured in a Vostok firn core (JM Barnola and L Arnaud, personal communication) are displayed in blue in Figure 5 alongside the "Lock-In" data in black. As explained in Section 3.2, these data were corrected for 5 cut bubbles using the "Lock-In" correction derived Section 3.1.2. In the case of Vostok, the Goujon et al. (2003) parametrization predicts closed porosity ratios (in red in Figure 5) that are more in line with the data than in the case of "Lock-In". Yet, an underestimation of the closed porosity subsists. However, it is possible that this discrepancy between the Goujon et al. (2003) parametrization and the Vostok data is due to a bad estimation of the cut-bubble effect. Comparison of the "Lock-In" and Vostok data suggests that the Vostok firn reaches closure at a higher density. This is corroborated by air content measurements: 10 the Vostok ice core has a lower total porous volume at isolation V i , also indicating that the firn reaches closure at a higher density. The difference between the two datasets is more pronounced than the difference between the "Lock-In" and Vostok Goujon et al. (2003) parametrizations in Figure 5. Therefore the temperature difference between the two sites does not fully explain the contrast between the two datasets, and it suggests that other climatic parameters such as the accumulation rate or the insolation also influences the porous network and the density range of bubble closure. This is consistent with the work of 15 Burr et al. (2018) that observed different pore network structures in the "Lock-In" and Dome C firns, two sites with similar temperatures.

Density as a predictor of a firn layer closure
In previous studies, the closed porosity ratio of a firn layer at a given polar site has primarily been linked to its density (Martinerie et al., 1992;Schwander et al., 1993;Mitchell et al., 2015;Burr et al., 2018), despite the firn heterogeneities. This 20 observation is confirmed by our work in the "Lock-In" firn, where both pycnometry and tomography indicate that the closed porosity ratio of a sample is primarily determined by its density. As depicted Figure S8 of the supplement, a larger dispersion is obtained when plotting the closed porosity ratios as a function of depth. In other words, two firn layers with same density, no matter how far apart they are in the firn column, have nearly the same closed porosity ratio. For example, a particular pycnometry sample taken at depth 97.40 m has a relative density of 0.90 and a closed porosity ratio close to 60%. Such values are 25 usually expected for samples around 102 m depth. Figure 5. Comparison of closed porosity ratio data and parametrizations. The "Lock-In" closed porosity ratio obtained by pycnometry and tomography are displayed as black hexagons (this study). The closed porosity ratios measured by pycnometry in the Vostok firn are displayed as blue circles (JM Barnola and L Arnaud, personal communication). The black curve is the Barnola parametrization for "Lock-In", using Vi = 0.114 cm 3 g −1 (derived from air content data). The orange and red curves are the "Lock-In" and Vostok parametrizations proposed by Goujon et al. (2003) using site temperatures as parameters. The green curve is the parametrization proposed for the B53 site by Schaller et al. (2017).
However, the fact that different firn layers exhibit a functional relationship between their density and closed porosity ratio, does not imply that they are equivalent in terms of pore structure. One can then wonder if especially high or low density layers display similar pore structure as the rest of the firn, or if they have a significantly different porous networks. To answer this question we used the tomography 3D images, which explicitly represent the porous network. Three different pore structure indicators were computed on each of the 1 cm thick slices.

5
The first indicator is the Structure Model Index or SMI. This index was developed by Hildebrand and Rüegsegger (1997) to analyze the structure of bones. It is traditionally used to quantify the degree of rodness, plateness, and spheriness in 3D images.
In this study the SMI will be used as a simple indicator of pore shape, considering that different pore morphologies should result in different SMI. As pointed out by Salmon et al. (2015) this index only provides sensible results when used on convex structures. For the "Lock-In" firn, the SMI was computed in each slice using the BoneJ plug-in within ImageJ (Doube et al., 10 2010). Results from the plug-in indicated that the porous structure is more than 90 % convex, justifying the use of SMI for firn studies (Gregory et al., 2014;Burr et al., 2018). Nonetheless, the obtained SMI values were sensitive to the resolution of the mesh used to compute the index.
The second indicator is the ratio of pore surface to pore volume, denoted S/V here. This surface to volume ratio was determined using the method presented in Krol and Löwe (2016), and is based on the determination of the two-point correlation function.
The S/V ratio is notably dependent on the size of the pore, as smaller objects tend to have higher surface to volume ratio.
Finally, we estimated the effective diffusivity of the firn samples (Schwander et al., 1988;Fabre et al., 2000;Freitag et al., 2002). The estimation was performed using the open-source matlab application TauFactor (Cooper et al., 2016). Effective diffusivity coefficients were obtained for vertical diffusion between the bottom and top sides of each 1 cm thick slice. However, 5 the resulting values should not be used to model gas diffusion in the firn, as the small size of the samples might not properly capture the large-scale behavior of the firn column (Fabre et al., 2000). Still, the comparison of effective diffusivity between individual samples can be used to highlight similarities or differences in their ability to transport gases. In that sense, our effective diffusivity coefficient is similar to the Qualitative Measure of Permeability (QMP) used by Fujita et al. (2009). Figure 6. Closed porosity ratio and morphology indicators against relative density: closed porosity ratio, S/V ratio, effective diffusivity, and structural model index. Each color corresponds to the depth where the firn samples were taken, and covers a 10 cm in the firn core.
In terms of general trend, the changes of the S/V ratio in Figure 6 are characterized by an increase with higher densities. This is to be expected as increasing S/V ratio is linked to a decrease in pore size with the compression of the firn. Moreover, our results are consistent with the ones reported by Burr et al. (2018) for the same "Lock-In" firn. They are also consistent with the S/V values reported by Gregory et al. (2014) for the WAIS Divide and Megadunes sites. Similarly, the decrease in effective diffusivity with density is consistent with the constriction of the pore network, where there is less and less straight 15 paths traversing the whole sample. Even though our results are not directly comparable to the ones of Freitag et al. (2002) due to the difference in methodology, we also find that the effective diffusivity is related to the open porosity by a power law.
However, the exponent found in our case is of 2.9, rather than the 2.1 reported by Freitag et al. (2002). Finally, the evolution of SMI is characterized by an increase in values with density, indicating a transformation of the pore shapes with increasing density. The values displayed in Figure 6 are consistent with the ones reported by Burr et al. (2018) for the "Lock-In" firn.
Similarly to Burr et al. (2018), our SMI values are higher than the ones reported by Gregory et al. (2014) for the WAIS Divide and Megadunes firns. This can be explained by the fact that the SkyScan software used by Gregory et al. (2014) yield lower SMI values (Salmon et al., 2015;Burr et al., 2018).

5
The tomography data of closed porosity ratio for the firn slices are also reported in Figure 6. Note that the closed porosity values displayed in Figure 6 have not been corrected for cut-bubbles. It is particularly apparent for the slices taken at the depth 99.85 m (in red Figure 6) that the firn stratification can lead to a wide range of closed porosity ratios within a 10 cm section.
Nonetheless, these layers with especially closed or open porosity do not appear as outliers in terms of pore morphology. Their morphological indicators are overall consistent with the other firn samples. Hence, this indicates that all layers tend to follow a 10 similar evolution with density not only in terms of closed porosity ratio, but also in terms of pore morphology and gas transport properties.

Origin of the density variability
We have seen in the previous sections that the degree of closure of a firn layer is primarily a function of its density. However, an important variability of the density is observed at the centimeter scale. One may then wonder about the origin of this density 15 variability in the trapping zone. Previous studies have advanced two mechanisms for the appearance of firn stratification in upper parts of the firn.
The first one is the softening of firn layers due to chemical impurities in the ice phase. Layers with a high concentration of specific ions densify faster, thus creating the density layering. Hörhold et al. (2012) observed a positive correlation between Ca 2+ and density anomalies in Greenland and Antarctic firn cores. They however do not argue that calcium is the ion specifi- hardening effect of NH + 4 , preventing the deformation and compaction of firn. The second stratification mechanism is detailed in Fujita et al. (2009), and takes into account the ice phase structure of firn layers. Briefly, due to surface metamorphism, some layers develop bonds between the ice clusters giving them more resistance to deformation. These layers, that are initially denser 25 at the surface, thus densify slower and are associated with negative density anomalies in deep firn. Using the high-resolution chemistry and density datasets of the "Lock-In" firn core we investigate the potential effect of chemistry on deep firn stratification, as proposed by Hörhold et al. (2012) and Fujita et al. (2016). The comparison of the high resolution density and liquid conductivity datasets indicates numerous patterns of covariation between the two signals, with denser layers having higher conductivities. Density and liquid conductivity variations for five one-meter long sections are displayed in Figure 7, with green and blue portions highlighting some of the observed covariations. This type of covariation is present throughout the whole trapping zone, and is consistent with the preferential deformation of firn with high ionic contents. Thus, it appears that the centimeter scale stratification observed in the trapping zone of "Lock-In" is in part related 5 to the presence of chemical impurities. Yet, liquid conductivity does not account for all the observed density variability. An illustration can be seen at the depth 90.30 m (pink section Figure 7), where a lower density layer does not show low values of liquid conductivity. This means one of two things. Either this less dense layer is due the absence of ice softening impurities but still contains enough conductive species for the liquid conductivity no to drop. In this case the density variability is still due to chemical impurities, but is simply not reflected in the liquid conductivity. Or this less dense layer is not due to a chemistry 10 based effect but rather to another mechanism, for example related to the ice structure as described by Fujita et al. (2009).
In order to better understand the effect of specific ions on the firn density, four one meter long sections have been measured using ion chromatography with a resolution of 2.5 cm. The concentration of H + was obtained by closing the ionic balance following Equation 3 of Legrand and Mayewski (1997). Since H + ions have a high mobility, H + concentration estimates could 15 help to explain the centimeter scale variations of liquid conductivity. The ion concentrations measured on the 113 to 114 m section are displayed in lower panels of Figure 8, with density and liquid conductivity variations above. On this section, the data highlight that in one occasion, the sulfate and conductivity undergo strong increases, uncorrelated with other species (green section in Figure 8). We interpret this spike as a volcanic event, depositing sulfuric acid. It is not associated with any specific effect on the firn density, and we conclude that sulfate has no impact on the layering in the trapping zone. Unfortunately, for 20 other species the results are inconclusive, as we were not able to decipher the effect of specific ions on densification. The same inconclusive observations were made on the three other sections (displayed in Figures S9 to S11 of the Supplement), where density variations could not be attributed to a particular species measured by ion chromatography. Measurements with a better spatial resolution could potentially help to decipher the impact of specific ion species, similarly to the work of Fujita et al.

Examples of stratified gas trapping
Thanks to the continuous methane record, we were able to identify layers which were not entirely closed below the bulk firnice transition. In the upper part of the methane record, these layers appear as spikes of high concentrations. This is due to  Figure 8. This still open layer is associated with low density and low liquid conductivity values. Please note that the methane axis is reversed in this figure, and in the following ones. The same observation was made on the ice core section displayed in Figure 9. Two late closure layers, depicted by the blue sections, are associated with lower density values and lower liquid conductivities than in the surrounding layers. All these examples are consistent with the 5 previous observation that less dense layers tend to be associated with low ionic contents and thus reach closure deeper in the firn. Figure 9 also reveals a low conductivity layer around 110.77m not associated with a low density value. This highlights that the presence of potential ice-softening ions does not systematically reflects in the liquid conductivity.
Similarly, an early closure layer was found in the bubbly-ice portion of the core ( Figure 10). Following Rhodes et al. (2016) and Fourteau et al. (2017), this layer was detected as a negative methane anomaly during a period of methane rise. This section 10 corresponds to a part of the core with a poor geometry due to a broken edge, and we were not able properly determine the absolute density values. Still, using the X-ray scan of the non broken part of the core we could identify the density variability.
However, as we do not know the actual thickness of the core that was crossed by the ray, the resulting density variability is expressed in arbitrary units and rather than in g cm −3 . The comparison between density and the methane records in Figure 10, indicates that the early closure layer corresponds to a denser layer. Air content measurements were performed on this section of the core. Results indicates that the early closure layer does not undergo outlier values of air content. This implies that the early closure layer trapped gases at similar densities as the rest of the firn. It is in line with the description of stratified gas trapping supported in this article, where dense layers close in advance but with the same critical density as others. Unfortunately, no chemical analysis could be performed on this section.

4 Conclusions
We studied the enclosure of gases at the East Antarctic site of "Lock-In", characterized by a temperature of −53.15 • C and an accumulation rate of 3.6 cm.we.yr −1 . The closed porosity ratio profile in the trapping zone was measured using two independent methods: pycnometry and tomography. Our findings suggest that the two methods yield consistent results when measuring similar samples, and can thus be interpreted in similar ways. However, pycnometry measurements are very sensitive 10 to the determination of sample volumes. On the other hand, tomography appears to provide more robust results and can be used to access more complex characterization of the firn, including pore size and shapes. One important remaining problem in the study of pore closure is the estimation of the cut-bubble effect in the trapping zone as it can drastically change the shape of the closed porosity curve and the estimation of the full closure density. Thus, we encourage the usage of large volume tomography for future studies of gas trapping and firn densification. Notably, large volume tomographies of firn where pycnometry data are Figure 10. Example of an early closure layer highlighted by the orange shade. The reduced methane concentration around 122.40 m is indicative of early gas trapping in the ice (reversed axis). This early closure layer is associated with high density values, but air content values similar with the adjacent layer. available could help understand the discrepancy between our results and the results of Schaller et al. (2017).
Using a gas trapping model, we found that under the assumption of iso-compression of open and closed pores, the measured air content and closed porosity ratios are inconsistent, with the model predicting to much enclosed air. The same overestimation of air content by the gas trapping model have been found in the Vostok ice core. We have shown that this discrepancy cannot be explained by a poor estimation of the closed porosity data or of the cut-bubble effect. On the other hand, the introduction of 5 a reduced compressibility for closed pores in the gas trapping model is able to resolve this discrepancy. Yet, this mechanism is not fully satisfactory as it requires a very limited compression of bubbles that has no direct supporting observations. The mechanism responsible for these low quantities of air trapped in ice cores is thus not clear.
Consistently with previous studies conducted at other sites, we observed a strong layering in the "Lock-In" firn, manifesting itself as centimeter scale density and firn structure variations (Hörhold et al., 2011;Fujita et al., 2016). Despite this heterogeneity, it appears that all firn layers roughly follow the same behavior relating their closed porosity ratio and pore structure to their density. Layers associated with positive density anomalies close in advance, but at the same critical density as others.
This vision simplifies the modeling of stratified gas trapping, as it reduces the firn to a stack of equivalent layers.
Finally, using liquid conductivity as a proxy for ionic content, our results suggest that when focusing on the trapping zone, 5 chemical impurities play a significant role in establishing the firn stratification. This does not mean that effects due to the ice phase structure, as for instance described by Fujita et al. (2009), does not strongly influence firn stratification in upper parts.
But, as a first approximation, the effect chemical impurities appears as the main mechanism for stratified gas trapping. Continuous impurity measurements along the trapping zone with a resolution higher than a centimeter could help decipher the species responsible for the presence of firn heterogeneities.

10
Code availability. The codes used in this study were developped using the Python3 language (with the readily available packages numpy, matplotlib, scipy, and skimage), and the macro language of ImageJ. The codes will be provided upon request to the corresponding authors.
Data availability. The high resolution density data in the trapping zone, the pycnometry and tomography closed porosity data, and the chemistry data (liquid conductivity and major ions) for Lock-In will be made accessible on the World Data Center for Paleoclimatology.
The Lock-In smoothed density and closed porosity data used in the gas trapping model will be made accessible as Supplementary material.

15
The Vostok pycnometry data will be published in a dedicated article. Meanwhile, they will be provided unpon request to the corresponding authors.
Author contributions. The high resolution density measurements were perfomed by CS, JF and KF. The chemistry measurements were carried out by RT, RM, LT and KF. The pycnometry system was restored by OM and LA, and the pycnometry measurements were performed by XF and KF. The tomography scans and segmentations were carried out by HL and MS. The air content measurements were performed by 20 VL. The Lock-In project was supervised by PM and XF. The code development for data processing and modeling was performed by KF and PM. All authors contributed to the discussion and interpretation of the data. The article was written by KF with the help of all co-authors. measurements, Jean-Robert Petit for his help with the solid conductivity data, Anouk Marsal for her help on the tomography data, Samuel Pion for his work on the pycnometry measurements, Alexis Burr for his help with the ImageJ macro language. We thank Émilie Capron for our useful discussions. We acknowledge the Automatic Weather Station Project run by Dr. Charles R. Stearns at the University of Wisconsin-Madison which is funded by the National Science Foundation of the United States of America for providing pressure data at Dome C. This work was supported by the French INSU/CNRS LEFE projects NEVE-CLIMAT and HEPIGANE.