Using ice core measurements from Taylor Glacier, Antarctica, to calibrate in situ cosmogenic 14C production rates by muons

Abstract

Cosmic rays entering the Earth's atmosphere produce showers of secondary particles such as protons, neutrons, and muons. The interaction of these particles with oxygen-16 (${}^{\mathrm{16}}\mathrm{O}$) in minerals such as ice and quartz can produce carbon-14 (${}^{\mathrm{14}}\mathrm{C}$). In glacial ice, ${}^{\mathrm{14}}\mathrm{C}$ is also incorporated through trapping of ${}^{\mathrm{14}}\mathrm{C}$-containing atmospheric gases (${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, ${}^{\mathrm{14}}\mathrm{CO}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$). Understanding the production rates of in situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ is important to deconvolve the in situ cosmogenic and atmospheric ${}^{\mathrm{14}}\mathrm{C}$ signals in ice, both of which contain valuable paleoenvironmental information. Unfortunately, the in situ ${}^{\mathrm{14}}\mathrm{C}$ production rates by muons (which are the dominant production mechanism at depths of $>\mathrm{6}$ m solid ice equivalent) are uncertain. In this study, we use measurements of in situ ${}^{\mathrm{14}}\mathrm{C}$ in ancient ice ($>\mathrm{50}$ ka) from the Taylor Glacier, an ablation site in Antarctica, in combination with a 2D ice flow model to better constrain the compound-specific rates of ${}^{\mathrm{14}}\mathrm{C}$ production by muons and the partitioning of in situ ${}^{\mathrm{14}}\mathrm{C}$ between CO${}_{\mathrm{2}}$, CO, and ${\mathrm{CH}}_{\mathrm{4}}$. Our measurements show that 33.7 % ($\pm \mathrm{11.4}\mathit{\%}$; 95 % confidence interval) of the produced cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ forms ${}^{\mathrm{14}}\mathrm{CO}$ and 66.1 % ($\pm \mathrm{11.5}\mathit{\%}$; 95 % confidence interval) of the produced cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ forms ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$. ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ represents a very small fraction ($<\mathrm{0.3}\mathit{\%}$) of the total. Assuming that the majority of in situ muogenic ${}^{\mathrm{14}}\mathrm{C}$ in ice forms ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, ${}^{\mathrm{14}}\mathrm{CO}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$, we also calculated muogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates that are lower by factors of 5.7 (3.6–13.9; 95 % confidence interval) and 3.7 (2.0–11.9; 95 % confidence interval) for negative muon capture and fast muon interactions, respectively, when compared to values determined in quartz from laboratory studies (Heisinger et al., 2002a, b) and in a natural setting (Lupker et al., 2015). This apparent discrepancy in muogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates in ice and quartz currently lacks a good explanation and requires further investigation.

1 Introduction

1.1

Potential applications of ${}^{\mathrm{14}}\mathrm{C}$ measurements in ice and in situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ production from ${}^{\mathrm{16}}\mathrm{O}$ in Earth's surface minerals

As snow accumulates on ice sheets, it gradually densifies into firn and ice (Herron and Langway, 1980). During the firn-to-ice transition, the air in the interstitial space between the ice grains becomes trapped into bubbles within the ice matrix (Buizert, 2013). Included in the paleoatmospheric air trapped in the bubbles are ${}^{\mathrm{14}}\mathrm{C}$-containing atmospheric gases (${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, ${}^{\mathrm{14}}\mathrm{CO}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$; Fireman and Norris, 1982). ${}^{\mathrm{14}}\mathrm{C}$ in ice is also produced through interactions of secondary cosmic rays with ${}^{\mathrm{16}}\mathrm{O}$ directly in the lattice of the ice grains (i.e., in situ; Lal et al., 1990). Following the cosmogenic nuclear reactions, the hot ${}^{\mathrm{14}}\mathrm{C}$ atom interacts with atoms in the surrounding ice lattice to produce ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, ${}^{\mathrm{14}}\mathrm{CO}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ (Lal et al., 1990; Petrenko et al., 2013).

Both the trapped atmospheric and in situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ signals in ice have unique applications. For example, the paleoatmospheric component of ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ in ice cores has been used to constrain past ${\mathrm{CH}}_{\mathrm{4}}$ emissions from old carbon reservoirs such as methane hydrates, permafrost, and geologic seeps (Dyonisius et al., 2020; Hmiel et al., 2020; Petrenko et al., 2009, 2017). Paleoatmospheric ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ can be potentially used for the absolute dating of ice core gases (Andree et al., 1984; Van De Wal et al., 1994) and to improve the radiocarbon calibration curve (Reimer et al., 2020; Hogg et al., 2020) in periods where tree ring data are not available. Measurements of ${}^{\mathrm{14}}\mathrm{CO}$ in the modern atmosphere have been used to constrain the oxidative capacity of the atmosphere (Brenninkmeijer et al., 1992; Petrenko et al., 2021), and thus, paleoatmospheric ${}^{\mathrm{14}}\mathrm{CO}$ in ice cores can be used for a similar application. The in situ cosmogenic component of ${}^{\mathrm{14}}\mathrm{CO}$ at ice core sites can be potentially be used to reconstruct the past cosmic ray flux (BenZvi et al., 2019). Finally, measurements of the in situ cosmogenic component of ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ and ${}^{\mathrm{14}}\mathrm{CO}$ can be used to constrain the accumulation/ablation rate of the ice core site (e.g., Lal et al., 1990; Lal and Jull, 1990). Unfortunately, the paleoatmospheric and in situ cosmogenic components of ${}^{\mathrm{14}}\mathrm{C}$ in ice exist in a combined form and cannot be separated analytically (Petrenko et al., 2016). To separate these signals, it is important to have accurate estimates of the cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates and the partitioning among the in situ produced ${}^{\mathrm{14}}\mathrm{C}$ species (${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, ${}^{\mathrm{14}}\mathrm{CO}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$) in ice.

In situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ production in ice is analogous to the production in quartz because both minerals share the same target atom (${}^{\mathrm{16}}\mathrm{O}$). Measurements of in situ cosmogenic nuclides (${}^{\mathrm{3}}\mathrm{He}$, ${}^{\mathrm{10}}\mathrm{Be}$, ${}^{\mathrm{14}}\mathrm{C}$, ${}^{\mathrm{21}}\mathrm{Ne}$, ${}^{\mathrm{26}}\mathrm{Al}$, and ${}^{\mathrm{36}}\mathrm{Cl}$) in near-surface rocks are commonly used as tools to constrain various Earth surface processes such as the timing of glacial retreat and erosion rates (Gosse and Phillips, 2001; Balco, 2020). Due to its short half-life of 5700 $\pm \mathrm{30}$ years (Kutschera, 2019), ${}^{\mathrm{14}}\mathrm{C}$ in quartz is uniquely suited to characterize surface processes on millennial timescales (e.g., Spector et al., 2019; Pendleton et al., 2019). In situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ measurements are also often paired with measurements of longer-lived nuclides such as ${}^{\mathrm{10}}\mathrm{Be}$ and ${}^{\mathrm{26}}\mathrm{Al}$ (e.g., Hippe, 2017; Skov et al., 2019) to study complex surface processes such as subglacial erosion and millennial-scale glacier retreats and/or re-advances.

In situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ in Earth's surface minerals is produced from ${}^{\mathrm{16}}\mathrm{O}$ by three nuclear reactions, namely (1) nucleon (neutron and proton)-induced spallation (Lal and Peters, 1967), (2) negative muon capture (Heisinger et al., 2002b), and (3) interactions with fast muons (Heisinger et al., 2002a). The depth dependence of the ${}^{\mathrm{14}}\mathrm{C}$ production rate for each mechanism in ice is shown in Fig. 1. Nucleon-induced spallation dominates the ${}^{\mathrm{14}}\mathrm{C}$ production at the surface but is quickly attenuated with depth. On the other hand, the relative contributions from the two muon mechanisms are lower near the surface but dominate at larger depths, as muons can penetrate deeper than nucleons (Fig. 1). Characterizing the in situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates from muons is especially important for applications of cosmogenic surface exposure dating where the samples might be exposed to the subsurface cosmic ray flux for an extended period. One example of this would be a bedrock that is covered by a relatively thin (e.g., tens of meters) glacier.

Figure 1

(a) In situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates scaled for the Taylor Glacier study site (77${}^{\circ }$44${}^{\prime }$ S, 162${}^{\circ }$10${}^{\prime }$ E; 526 m elevation) from the three nuclear mechanisms, i.e., nucleon-induced spallation ($P_{\mathrm{n}}$), negative muon capture ($P_{\mathrm{neg}}$), and fast muon interactions ($P_{\mathrm{fast}}$). (b) Depth profiles of the three production mechanisms normalized to their respective surface production rates (the respective surface production rates are shown in the legend). For ${}^{\mathrm{14}}\mathrm{C}$ production from neutron spallation, we used the surface production rate estimate from Young et al. (2014) with scaling from the Lifton et al. (2014) LSDn model. For the two muon mechanisms (negative and fast muons), we used the production rate model from Balco et al. (2008), which follows parameterizations by Heisinger et al. (2002a, b).

[Figure omitted. See PDF]

Understanding the muogenic ${}^{\mathrm{14}}\mathrm{C}$ component is especially important for ${}^{\mathrm{14}}\mathrm{C}$ studies in ice. Prior studies have shown that, at snow accumulation sites, most of the in situ ${}^{\mathrm{14}}\mathrm{C}$ produced in the firn (including the majority of neutron-produced ${}^{\mathrm{14}}\mathrm{C}$) is lost to the atmosphere via gas movement in the firn open porosity (Petrenko et al., 2013; Van der Kemp et al., 2000; Wilson and Donahue, 1990). In situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ mainly starts to accumulate in deeper ice where gas exchange with the atmosphere no longer happens, and at these depths, the ${}^{\mathrm{14}}\mathrm{C}$ production is entirely from the muon mechanisms. Thus, the in situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ signal in traditional deep ice cores is dominated by production from muons, and constraining the muogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates is critical to disentangle the in situ cosmogenic and atmospheric ${}^{\mathrm{14}}\mathrm{C}$ signals in ice cores. Unfortunately, the in situ ${}^{\mathrm{14}}\mathrm{C}$ production rates by muons in both ice and quartz are still highly uncertain (Hippe, 2017).

The production rates of cosmogenic nuclides are usually determined from calibration sites where independent controls on exposure history are available, such as ${}^{\mathrm{14}}\mathrm{C}$ dating from organic materials (e.g., Lifton et al., 2015) or argon (${}^{\mathrm{40}}\mathrm{Ar}{/}^{\mathrm{39}}\mathrm{Ar}$) dating from lava flows (e.g., Balbas and Farley, 2020; Fenton et al., 2019). However, the commonly used estimates of muogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates (for both negative muon capture and fast muon reactions) were derived through the laboratory irradiation of artificial target compounds (Heisinger et al., 2002a, b). To our knowledge, there is only one prior study (Lupker et al., 2015) that provided estimates of total muogenic in situ ${}^{\mathrm{14}}\mathrm{C}$ production rates based on measurements in a natural setting. Using ${}^{\mathrm{14}}\mathrm{C}$ measurements from a 15.5 m deep quartzite core from Leymon High, Spain, Lupker et al. (2015) estimated a sea level high latitude (SLHL) surface production rate of 3.34 ($+\mathrm{0.43}/-\mathrm{1.07}$) ${}^{\mathrm{14}}\mathrm{C}$ atoms per gram of quartz per year for negative muon capture and 0 ($+\mathrm{0.42}/-\mathrm{0.00}$) ${}^{\mathrm{14}}\mathrm{C}$ atoms per gram of quartz per year for fast muon interactions ($\mathrm{1}\mathit{\sigma }$ uncertainties). The large uncertainties in the ${}^{\mathrm{14}}\mathrm{C}$ production rates (especially the production rate from fast muons) estimated by Lupker et al. (2015) were due to relatively large measurement uncertainty for their deepest samples and a small contribution to the ${}^{\mathrm{14}}\mathrm{C}$ signal from fast muons. Balco (2017) also refitted the Leymon High data and obtained similar results regarding the magnitude of the SLHL ${}^{\mathrm{14}}\mathrm{C}$ production rate from negative muon capture.

Petrenko et al. (2016) used ${}^{\mathrm{14}}\mathrm{C}$ measurements (${}^{\mathrm{14}}\mathrm{CO}$, ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$) in $>\mathrm{50}$ ka ice for the 2–20 m depth range from Taylor Glacier, Antarctica, to constrain the ${}^{\mathrm{14}}\mathrm{C}$ production rates in ice. The old age of the ice ensured that all in situ cosmogenic and paleoatmospheric ${}^{\mathrm{14}}\mathrm{C}$ inherited from the ice accumulation site had decayed away. Unfortunately, Petrenko et al. (2016) were unable to accurately constrain the total ${}^{\mathrm{14}}\mathrm{C}$ production rates because of the high uncertainty resulting from the melt–extraction technique used to obtain their ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements (see Sect. 1.3).

1.2

Overview of ${}^{\mathrm{14}}\mathrm{C}$ production from muons

Following Heisinger et al. (2002b), the production rate of ${}^{\mathrm{14}}\mathrm{C}$ (atoms per gram per year) by negative muon capture ($P_{\mathrm{neg}}$) as a function of lithospheric depth ($h$; typically in g cm${}^{-\mathrm{2}}$) is given by the following: $$\begin{array}{}\text{1}& {\displaystyle }{P}_{\mathrm{neg}}\left(h\right)=R_{\mathit{\mu }-}\left(h\right)\cdot f_{\mathrm{tot}}\text{2}& {\displaystyle }{f}_{\mathrm{tot}}=f_{\mathrm{C}}\cdot f_{\mathrm{D}}\cdot f^{*},\end{array}$$ where $R_{\mathit{\mu }-}\left(z\right)$ is the stopping rate of negative muons (muons per gram per year) at lithospheric depth $h$, and $f_{\mathrm{tot}}$ is the overall probability of ${}^{\mathrm{14}}\mathrm{C}$ production in ice from a stopped negative muon (unitless). The stopping rate of negative muons at the given depth $R_{\mathit{\mu }-}\left(h\right)$ has been empirically determined from measurements at deep underground laboratories (Heisinger et al., 2002b). The lithospheric depth ($h$) is a product of actual depth ($z$) and density ($\mathit{\rho }$) of the target mineral (${\mathit{\rho }}_{\mathrm{ice}}=\mathrm{0.92}$ g cm${}^{-\mathrm{3}}$).

The total probability ($f_{\mathrm{tot}}$) of ${}^{\mathrm{14}}\mathrm{C}$ production from negative muon capture is expressed by the product of the chemical compound factor ($f_{\mathrm{C}}$), representing the probability that the stopped muon is captured by one of the target atoms (${}^{\mathrm{16}}\mathrm{O}$ in case of ${}^{\mathrm{14}}\mathrm{C}$ production), the probability that the negative muon does not decay in the K shell before nuclear capture ($f_{\mathrm{D}}$), and the effective probability for the production of cosmogenic nuclides after $\mathit{\mu }$- capture by the target atom ($f^{*}$; Eq. 2; Heisinger et al., 2002b; Lupker et al., 2015). All probability ($f$) terms in Eq. (2) are unitless. From experiments involving laboratory irradiation of artificial targets, the overall probability ($f_{\mathrm{tot}}$) for ${}^{\mathrm{14}}\mathrm{C}$ production in ice from negative muon was estimated to be $\mathrm{0.025}\pm \mathrm{0.002}$ (Heisinger et al., 2002b).

An expression for the production rate of nuclides by fast muon interactions ($P_{\mathrm{fast}}$) as a function of lithospheric depth ($h$) is given by Heisinger et al. (2002a): $$\begin{array}{}\text{3}& {\displaystyle }{P}_{\mathrm{fast}}\left(h\right)={\mathit{\sigma }}_{\mathrm{0}}\cdot \mathit{\beta }\left(h\right)\mathit{\varphi }\left(h\right)\cdot \stackrel{\mathrm{‾}}{E}\left(h\right)\mathit{\alpha }\cdot N\text{4}& {\displaystyle }\mathit{\beta }\left(h\right)=\mathrm{0.846}-\mathrm{0.015}\ln (h+\mathrm{1})+\mathrm{0.003139}\left(\ln \right(h+\mathrm{1}){)}^{\mathrm{2}},\end{array}$$ where $\mathit{\varphi }\left(h\right)$ is the total muon flux at depth $z$ (muons per centimeter squared per year per steradian), ${\mathit{\sigma }}_{\mathrm{0}}$ is the reference nuclear reaction cross section at muon energy of 1 GeV (millibarn, mb), $\mathit{\beta }\left(h\right)$ is the unitless parameterized depth-dependence factor (Eq. 4), $\stackrel{\mathrm{‾}}{E}\left(h\right)$ is the mean muon energy at depth $h$ (GeV), $\mathit{\alpha }$ is a power factor that describes the energy dependence of the cross section (unitless; $\mathit{\alpha }=\mathrm{0.75}$), and $N$ is the number of target nuclei per gram of the target mineral. The overall production rate of ${}^{\mathrm{14}}\mathrm{C}$ from fast muons provided by Heisinger et al. (2002a) has a high ($\pm \mathrm{50}\mathit{\%}$) uncertainty because of the uncertainty in the reference nuclear reaction cross section ${\mathit{\sigma }}_{\mathrm{0}}$ (${\mathit{\sigma }}_{\mathrm{0}}=\mathrm{0.0088}\pm \mathrm{0.0049}$ mb). Following Lupker et al. (2015), in this study we used $f_{\mathrm{tot}}$ and ${\mathit{\sigma }}_{\mathrm{0}}$ as tuning parameters for the two muogenic production mechanisms in a cosmogenic nuclide production model (Sect. 3.2) to fit our ${}^{\mathrm{14}}\mathrm{C}$ measurements.

1.3

Gas extraction methods for ice core ${}^{\mathrm{14}}\mathrm{C}$ analysis

Common methods to liberate gas trapped in ice core bubbles include melting (wet extraction; e.g., Sowers et al., 1992; Mitchell et al., 2011) and mechanical destruction of the ice lattice (dry extraction; e.g., Bereiter et al., 2013; Ahn et al., 2009; Zumbrunn et al., 1982). Dry extraction is generally preferable for CO${}_{\mathrm{2}}$ analysis because the presence of liquid water in a wet extraction introduces extraneous CO${}_{\mathrm{2}}$ from the carbonate–acid reaction between the meltwater and impurities in the ice (e.g., Delmas et al., 1980; Raynaud et al., 1982). Multiple studies of ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ in ice have used dry extraction methods (e.g., Van De Wal et al., 1994, 2007; Smith et al., 2000; Van der Kemp et al., 2000). However, dry extraction systems (e.g., Lüthi et al., 2008) can potentially introduce biases in the CO${}_{\mathrm{2}}$ mole fraction [CO${}_{\mathrm{2}}$] due to incomplete gas extraction (Bereiter et al., 2015). Considering that the in situ cosmogenic production of ${}^{\mathrm{14}}\mathrm{C}$ occurs directly in the ice lattice (Lal et al., 1990), it has been argued that dry extraction may also not liberate all of the ${}^{\mathrm{14}}\mathrm{C}$ from the ice (e.g., van Roijen et al., 1994).

Other studies of ${}^{\mathrm{14}}\mathrm{C}$ in ice (e.g., Lal et al., 1990, 1997, 2001; Jull et al., 1994) have used wet extraction methods. These wet extraction studies involved an addition of acid to drive off all dissolved CO${}_{\mathrm{2}}$ from the meltwater (Lal et al., 1990, 1997, 2001; Jull et al., 1994). The acidification process may have resulted in an additional CO${}_{\mathrm{2}}$ release from impurities in the ice (e.g., carbonate dust). In dust-rich Greenland ice, the presence of liquid water in a wet extraction produced excess ${\mathrm{CH}}_{\mathrm{4}}$ artifacts (Lee et al., 2020). It is thus possible that a wet extraction approach for ${}^{\mathrm{14}}\mathrm{C}$ analysis may also result in additional C release from organics in the ice which are not ${}^{\mathrm{14}}\mathrm{C}$ free.

A third method to liberate gases trapped in ice cores is sublimation under vacuum (e.g., Wilson and Donahue, 1989, 1990; Wilson and Long, 1997; Siegenthaler et al., 2005; Schmitt et al., 2011). Sublimation can occur when the pressure and temperature on the surface of the ice are below the triple point of the water-phase-change diagram. In addition to being free of problems associated with wet extraction methods, sublimation guarantees 100 % gas extraction efficiency (Schmitt et al., 2011; Bereiter et al., 2013, 2015), which includes any ${}^{\mathrm{14}}\mathrm{C}$ trapped in the ice lattice. Therefore, sublimation is likely an optimal method for ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements in ice.

This study presents new ${}^{\mathrm{14}}\mathrm{C}$ measurements in three gas species (${}^{\mathrm{14}}\mathrm{CO}$, ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$) in ancient ($>\mathrm{50}$ ka) ice from the ablation zone of the Taylor Glacier, Antarctica, to constrain the compound-specific ${}^{\mathrm{14}}\mathrm{C}$ production rates in ice by muons. Ice at this location does not contain a significant amount of ${}^{\mathrm{14}}\mathrm{C}$ inherited from the accumulation site (Petrenko et al., 2016), and the ${}^{\mathrm{14}}\mathrm{C}$ content is entirely due to production by muons during transport within the glacier. We improved on the earlier work by Petrenko et al. (2016) by (1) using a newly developed ice sublimation extraction device for ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements (see Sect. 2.3.2), (2) collecting deeper samples to $\sim \mathrm{72}$ m to better characterize the ${}^{\mathrm{14}}\mathrm{C}$ production rate from the fast muon mechanism, and (3) using a more realistic 2D ice flow model from Buizert et al. (2012) to account for the flow trajectory and exposure history of the samples (see Sect. 3.1).

2 Field sampling and analytical methods

2.1 Site description

The blue ice area of the Taylor Glacier (Fig. 2) provides access to near-unlimited amounts of well-dated ancient ice (Baggenstos et al., 2017; Bauska et al., 2016; Menking et al., 2019; Schilt et al., 2014; Shackleton et al., 2020). This allows ice from Taylor Glacier to be measured for ultra-trace-gas species that require a very large amount of ice (Dyonisius et al., 2020; Petrenko et al., 2016, 2017; Buizert et al., 2014). In this study, we used the same site as Petrenko et al. (2016; 77${}^{\circ }$43.699${}^{\prime }$ S, 161${}^{\circ }$43.179${}^{\prime }$), where ice $>\mathrm{50}$ ka in age at the surface has been previously identified.

2.2 Field sampling

Approximately 1000 kg of ice is needed to obtain both the necessary ${\mathrm{CH}}_{\mathrm{4}}$-derived and CO-derived C mass for ${}^{\mathrm{14}}\mathrm{C}$ analyses. Because of this large sample requirement, and to avoid post-coring in situ ${}^{\mathrm{14}}\mathrm{C}$ production at the surface, the melt extraction for ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ and ${}^{\mathrm{14}}\mathrm{CO}$ samples was performed on site using the large-volume melter apparatus and technique described in Petrenko et al. (2016). The liberated air was transferred to 34.9 L electropolished stainless steel canisters and shipped to our laboratories for processing and analyses. Similar to other studies using this large-volume ice melter (e.g., Dyonisius et al., 2020; Petrenko et al., 2016, 2017), four procedural blanks (two with modern ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ standard gas and two with ${}^{\mathrm{14}}\mathrm{C}$ dead ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ standard gas) were collected in the field. These field procedural blanks allow us to characterize the addition of extraneous ${}^{\mathrm{14}}\mathrm{C}$ to the samples. The standard gases used in the field procedural blanks were passed through a Sofnocat 423 reagent, which removes CO (and thus ${}^{\mathrm{14}}\mathrm{CO}$) but leaves ${\mathrm{CH}}_{\mathrm{4}}$ (and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$) intact.

Figure 2

Map of the Taylor Glacier study site. The sampling location is marked by a red cross on the map. The orange star in the inset map shows the location of the Taylor Glacier relative to the Antarctic continent. The map was made using ArcGIS Pro, with imagery layers from Esri and Earthstar Geographics.

[Figure omitted. See PDF]

The sampling scheme for this study is shown in Fig. S1 in the Supplement. We used the 9.5 in. (241 mm) diameter Blue Ice Drill (BID; Kuhl et al., 2014) to collect seven large-volume samples during the 2015/2016 austral summer field season for ${}^{\mathrm{14}}\mathrm{CO}$ and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ analyses. The surface sample was collected from $\mathrm{21}\times \mathrm{1.5}$ m deep shallow cores, each with an average mid-depth of $\sim \mathrm{0.75}$ m. Six additional deep samples with mid-depths of 19.5, 30, 40.5, 51, 61.5, and 72 m were also collected by combining ice from three $\sim \mathrm{78}$ m deep boreholes. Each of the deep large-volume samples spanned approximately 10.5 m depth. Continuous sticks of ice subsamples ($\mathrm{3}\times \mathrm{3}$ cm; spanning the whole length of the core) were taken from one of the three ice core boreholes (TGDeep3) for age control (see Sect. S3 in the Supplement). The continuous sample sticks were measured for the ${\mathrm{CH}}_{\mathrm{4}}$ mole fraction [${\mathrm{CH}}_{\mathrm{4}}$] using the continuous flow analysis (CFA) system described in Rhodes et al. (2013) at Oregon State University (OSU).

In addition to the large-volume samples, we collected 26 smaller subsamples ($\sim \mathrm{1.5}$–2 kg) from 13 depth levels and two boreholes for ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements. Each depth level contained a pair of replicates; however, only 9 out of the 13 replicate pairs were true replicates (i.e., collected from the same borehole and cut from the same depth interval; Fig. S1). Collecting same depth-adjacent samples below 50 m depth from a single borehole was challenging because of the reduced core quality (i.e., more fractures in the ice), and thus, the replicates had to be collected from a different borehole. Immediately after removal from the borehole, ice samples become exposed to a more intense cosmic ray bombardment (post-coring in situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ production). Five artificial bubble-free ice (BFI) samples were manufactured in the field, following methods from Mitchell et al. (2011), but upscaled to produce 1.5–2 kg samples. The field-produced BFI samples were shipped together with the collected glacial ice samples to characterize the effects of the post-coring in situ cosmogenic ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ production in the samples.

2.3 Laboratory analytical methods

2.3.1

Large-volume samples for ${}^{\mathrm{14}}\mathrm{CO}$ and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ measurements

The detailed approach for sample processing, measurements, and associated procedural corrections for the large-volume samples have been previously described in detail (Petrenko et al., 2016). In this section, we only provide a brief overview and highlight the differences between our methods and those of Dyonisius et al. (2020). First, the ${\mathit{\delta }}^{\mathrm{13}}$${\mathrm{CH}}_{\mathrm{4}}$ measurements were conducted at the Institute of Arctic and Alpine Research (INSTAAR), following the methods described by Miller et al. (2002; Table S1 in the Supplement). The ${\mathit{\delta }}^{\mathrm{13}}$${\mathrm{CH}}_{\mathrm{4}}$ measurements were not corrected for gravitational (Sowers et al., 1992) and diffusive isotopic fractionation (Buizert et al., 2013) because these corrections are only necessary to reconstruct the paleoatmospheric ${\mathit{\delta }}^{\mathrm{13}}$${\mathrm{CH}}_{\mathrm{4}}$ signal. In this study, the ${\mathit{\delta }}^{\mathrm{13}}$${\mathrm{CH}}_{\mathrm{4}}$ values are only used to normalize and calculate the absolute ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ abundance (in molec. per gram of ice).

The large-volume samples and field procedural blanks were measured for [${\mathrm{CH}}_{\mathrm{4}}$] using a gas chromatograph–multidetector (GC-MD) system (Prinn et al., 2008; Table S2). Pressure in the sample canisters was measured using a Paroscientific, Inc., Digiquartz^® series 740 absolute pressure transducer at the Scripps Institution of Oceanography (SIO) for total air content (TAC) determination (Table S3). Two of the field procedural blanks were also measured for the $\mathrm{Kr}/{\mathrm{N}}_{\mathrm{2}}$, $\mathrm{Xe}/{\mathrm{N}}_{\mathrm{2}}$, and $\mathrm{Xe}/\mathrm{Kr}$ ratio (Table S4) at the Scripps Institution of Oceanography (SIO), following the procedures described in Bereiter et al. (2018). The noble gas ratios were used to constrain the degree of gas solubility during the melt extraction. The large-volume samples were measured for CO mole fraction [CO] using a Picarro G2401 analyzer (Table S5) and again for pressure at the University of Rochester (UR; Table S4).

The ${\mathrm{CH}}_{\mathrm{4}}$ in the large-volume samples and blanks was combusted to CO${}_{\mathrm{2}}$, cryogenically separated, and flame-sealed in glass ampules using the air processing line at the University of Rochester (Dyonisius et al., 2020). We also processed $\mathrm{3}\times \mathrm{100}$ $\mathrm{\mu }$g of ${\mathrm{CH}}_{\mathrm{4}}$-derived C samples each from the modern ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ standard gas and ${}^{\mathrm{14}}\mathrm{C}$ dead standard gas used for the field procedural blanks. The sample air that remained after ${\mathrm{CH}}_{\mathrm{4}}$ processing ($\sim \mathrm{10}$ L standard temperature pressure, STP) was diluted with a gas containing $\mathrm{10.02}\pm \mathrm{0.26}$ $\mathrm{\mu }$mol mol${}^{-\mathrm{1}}$ (95 % confidence interval, CI) of ${}^{\mathrm{14}}\mathrm{C}$-depleted CO (${}^{\mathrm{14}}\mathrm{CO}$ $=\mathrm{0.19}\pm \mathrm{0.08}$ pMC; 95 % CI) to increase the CO-derived C mass for the accelerator mass spectrometry (AMS) measurements. The dilutant gas was measured for ${\mathit{\delta }}^{\mathrm{13}}$CO using methods described in Vimont (2017; ${\mathit{\delta }}^{\mathrm{13}}$CO $=-\mathrm{23.36}\pm \mathrm{0.2}$ ‰; 95 %CI).

The CO- and ${\mathrm{CH}}_{\mathrm{4}}$-derived CO${}_{\mathrm{2}}$ was graphitized using the Australian Nuclear Science and Technology Organization (ANSTO) micro-furnaces, following Yang and Smith (2017). We used the ${}^{\mathrm{14}}\mathrm{C}$ activity measured on the 100 $\mathrm{\mu }$g samples as the true ${}^{\mathrm{14}}\mathrm{C}$ activity of the standard gases (Table S6). Because of the larger sample size, the effect of extraneous C introduced by graphitization on these 100 $\mathrm{\mu }$g samples is assumed to be negligible. Using a mass balance approach described in Petrenko et al. (2017), the total extraneous C mass for the ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ samples was determined to be $\mathrm{0.63}\pm \mathrm{0.28}$ $\mathrm{\mu }$gC, and the corresponding ${}^{\mathrm{14}}\mathrm{C}$ activity for the extraneous C was $\mathrm{16.7}\pm \mathrm{10.2}$ pMC (95 % CI).

Table 1

Measured ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, ${}^{\mathrm{14}}\mathrm{CO}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ after all associated corrections and calculated total ${}^{\mathrm{14}}\mathrm{C}$, ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}{/}^{\mathrm{14}}\mathrm{CO}$ ratios, ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, and ${}^{\mathrm{14}}\mathrm{CO}$ fractions. All errors presented indicate the 95 % CI. NA – not available.

Mid-depth	${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$	${}^{\mathrm{14}}\mathrm{CO}$	${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$	Total ${}^{\mathrm{14}}\mathrm{C}$	${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}{/}^{\mathrm{14}}\mathrm{CO}$	${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$	${}^{\mathrm{14}}\mathrm{CO}$
(m)	(molec. per gram	(molec. per gram	(molec. per gram	(atoms per gram	ratio	fraction	fraction
	of ice)	of ice)	of ice)	of ice)
2.25	$\mathrm{145.5}\pm \mathrm{32.0}$	$\mathrm{45.2}\pm {\mathrm{3.2}}^{\mathrm{a}}$	$\mathrm{0.418}\pm {\mathrm{0.052}}^{\mathrm{a}}$	$\mathrm{191.1}\pm \mathrm{32.2}$	$\mathrm{0.0092}\pm \mathrm{0.0013}$	$\mathrm{0.76}\pm \mathrm{0.18}$	$\mathrm{0.24}\pm \mathrm{0.04}$
	$\mathrm{123.3}\pm \mathrm{28.5}$			$\mathrm{168.9}\pm \mathrm{28.7}$		$\mathrm{0.73}\pm \mathrm{0.18}$	$\mathrm{0.27}\pm \mathrm{0.05}$
3.65	$\mathrm{88.5}\pm \mathrm{7.8}$	$\mathrm{41.9}\pm {\mathrm{2.2}}^{\mathrm{a}}$	$\mathrm{0.327}\pm {\mathrm{0.026}}^{\mathrm{a}}$	$\mathrm{130.7}\pm \mathrm{8.1}$	$\mathrm{0.0078}\pm \mathrm{0.0007}$	$\mathrm{0.68}\pm \mathrm{0.07}$	$\mathrm{0.32}\pm \mathrm{0.03}$
	$\mathrm{98.0}\pm \mathrm{7.8}$			$\mathrm{140.2}\pm \mathrm{8.1}$		$\mathrm{0.70}\pm \mathrm{0.07}$	$\mathrm{0.30}\pm \mathrm{0.02}$
6.85	$\mathrm{64.2}\pm \mathrm{7.8}$	$\mathrm{36.4}\pm {\mathrm{1.8}}^{\mathrm{a}}$	$\mathrm{0.273}\pm {\mathrm{0.020}}^{\mathrm{a}}$	$\mathrm{100.9}\pm \mathrm{8.0}$	$\mathrm{0.0075}\pm \mathrm{0.0007}$	$\mathrm{0.64}\pm \mathrm{0.08}$	$\mathrm{0.36}\pm \mathrm{0.03}$
	$\mathrm{69.4}\pm \mathrm{7.8}$			$\mathrm{106.1}\pm \mathrm{8.0}$		$\mathrm{0.65}\pm \mathrm{0.08}$	$\mathrm{0.34}\pm \mathrm{0.03}$
10${}^{\mathrm{b}}$	$\mathrm{50.6}\pm \mathrm{7.8}$	$\mathrm{31.4}\pm {\mathrm{1.6}}^{\mathrm{a}}$	NA	$\mathrm{82.2}\pm {\mathrm{8.0}}^{\mathrm{b}}$	NA	$\mathrm{0.62}\pm \mathrm{0.10}$	$\mathrm{0.38}\pm \mathrm{0.04}$
	$\mathrm{54.3}\pm \mathrm{7.8}$			$\mathrm{86}\pm {\mathrm{7.9}}^{\mathrm{b}}$		$\mathrm{0.63}\pm \mathrm{0.10}$	$\mathrm{0.37}\pm \mathrm{0.04}$
15	$\mathrm{60.9}\pm \mathrm{7.8}$	$\mathrm{26.9}\pm {\mathrm{1.4}}^{\mathrm{a}}$	$\mathrm{0.206}\pm {\mathrm{0.016}}^{\mathrm{a}}$	$\mathrm{88.0}\pm \mathrm{7.9}$	$\mathrm{0.0077}\pm \mathrm{0.0007}$	$\mathrm{0.69}\pm \mathrm{0.10}$	$\mathrm{0.31}\pm \mathrm{0.03}$
	$\mathrm{54.6}\pm \mathrm{7.8}$			$\mathrm{81.7}\pm \mathrm{7.9}$		$\mathrm{0.67}\pm \mathrm{0.10}$	$\mathrm{0.33}\pm \mathrm{0.04}$
19.5	$\mathrm{52.4}\pm \mathrm{7.8}$	$\mathrm{23.9}\pm {\mathrm{1.2}}^{\mathrm{a}}$	$\mathrm{0.182}\pm {\mathrm{0.016}}^{\mathrm{a}}$	$\mathrm{76.5}\pm \mathrm{7.9}$	$\mathrm{0.0076}\pm \mathrm{0.0008}$	$\mathrm{0.69}\pm \mathrm{0.11}$	$\mathrm{0.31}\pm \mathrm{0.04}$
	$\mathrm{49.6}\pm \mathrm{7.8}$			$\mathrm{73.7}\pm \mathrm{7.9}$		$\mathrm{0.67}\pm \mathrm{0.11}$	$\mathrm{0.32}\pm \mathrm{0.04}$
40.5	$\mathrm{36.4}\pm \mathrm{7.8}$	$\mathrm{15.8}\pm \mathrm{1.0}$	$\mathrm{0.119}\pm \mathrm{0.013}$	$\mathrm{52.3}\pm \mathrm{7.8}$	$\mathrm{0.0075}\pm \mathrm{0.0010}$	$\mathrm{0.70}\pm \mathrm{0.15}$	$\mathrm{0.30}\pm \mathrm{0.05}$
	$\mathrm{37.2}\pm \mathrm{7.8}$			$\mathrm{53.1}\pm \mathrm{7.9}$		$\mathrm{0.70}\pm \mathrm{0.15}$	$\mathrm{0.30}\pm \mathrm{0.05}$
51	$\mathrm{31.1}\pm \mathrm{7.8}$	$\mathrm{13.2}\pm \mathrm{0.9}$	$\mathrm{0.097}\pm \mathrm{0.014}$	$\mathrm{44.4}\pm \mathrm{7.8}$	$\mathrm{0.0073}\pm \mathrm{0.0012}$	$\mathrm{0.70}\pm \mathrm{0.18}$	$\mathrm{0.30}\pm \mathrm{0.06}$
	$\mathrm{28.0}\pm \mathrm{7.8}$			$\mathrm{41.3}\pm \mathrm{7.8}$		$\mathrm{0.68}\pm \mathrm{0.19}$	$\mathrm{0.32}\pm \mathrm{0.06}$
61.5	$\mathrm{22.1}\pm \mathrm{7.8}$	$\mathrm{11.3}\pm \mathrm{0.7}$	$\mathrm{0.079}\pm \mathrm{0.013}$	$\mathrm{33.5}\pm \mathrm{7.8}$	$\mathrm{0.0070}\pm \mathrm{0.0012}$	$\mathrm{0.66}\pm \mathrm{0.24}$	$\mathrm{0.34}\pm \mathrm{0.08}$
	$\mathrm{31.0}\pm \mathrm{7.8}$			$\mathrm{42.4}\pm \mathrm{7.8}$		$\mathrm{0.73}\pm \mathrm{0.19}$	$\mathrm{0.27}\pm \mathrm{0.05}$
72	$\mathrm{11.8}\pm \mathrm{7.8}$	$\mathrm{10.8}\pm \mathrm{0.7}$	$\mathrm{0.080}\pm \mathrm{0.013}$	$\mathrm{22.7}\pm \mathrm{7.8}$	$\mathrm{0.0074}\pm \mathrm{0.0013}$	$\mathrm{0.52}\pm \mathrm{0.34}$	$\mathrm{0.48}\pm \mathrm{0.17}$
	$\mathrm{13.6}\pm \mathrm{7.8}$			$\mathrm{24.5}\pm \mathrm{7.8}$		$\mathrm{0.55}\pm \mathrm{0.32}$	$\mathrm{0.44}\pm \mathrm{0.14}$

${}^{\mathrm{a}}$ Data from Petrenko et al. (2016). ${}^{\mathrm{b}}$ The total ${}^{\mathrm{14}}\mathrm{C}$ value for the 10 m sample was scaled by a factor of $\mathrm{1.003}\pm \mathrm{0.003}$ (95 % CI) to account for the lack of ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ measurements (Sect. 4.2).

Figure 3

(a) Measured ${}^{\mathrm{14}}\mathrm{CO}$ molec. per gram of ice after all corrections. (b) Measured ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ molec. per gram of ice after all corrections. (c) Measured ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ molec. per gram of ice after all corrections. (d) Total measured ${}^{\mathrm{14}}\mathrm{C}$ atoms per gram of ice. This represents the sum of ${}^{\mathrm{14}}\mathrm{CO}$, ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$, and sublimation-based ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements. All error bars shown in this figure are 95 % CI.

[Figure omitted. See PDF]

In prior studies (e.g., Dyonisius et al., 2020; Petrenko et al., 2017), ${}^{\mathrm{14}}\mathrm{CO}$ measurements from the field procedural blanks were used to characterize the effects of extraneous ${}^{\mathrm{14}}\mathrm{C}$ addition from sample extraction, handling, storage, transport, and processing (including the graphitization step). For this study, the field procedural blanks were still used to characterize the effects from in situ production of ${}^{\mathrm{14}}\mathrm{CO}$ in the sample air canisters by cosmic rays during storage and transport. However, to better characterize the effects from the addition of extraneous C during the graphitization process, we used a linear empirical correction from 10 commensurately sized ${}^{\mathrm{14}}\mathrm{C}$ standards and blanks at ANSTO (see Fig. S2a; Table S7), following Petrenko et al. (2021). This approach has the benefit of bracketing the effects of extraneous C from graphitization at ANSTO with low and high ${}^{\mathrm{14}}\mathrm{C}$ standards, similar to the approach for the ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ samples. The ${}^{\mathrm{14}}\mathrm{CO}$ blank for this sample set is $\mathrm{22.45}\pm \mathrm{3.24}$ molec. ${}^{\mathrm{14}}\mathrm{CO}$ per cc of air STP (95 % CI), which is higher than the ${}^{\mathrm{14}}\mathrm{CO}$ blanks reported in Dyonisius et al. (2020). This is mainly because there was an extra year between the retrieval and processing of the samples (thus, there was more in situ ${}^{\mathrm{14}}\mathrm{CO}$ production in sample canisters during storage). The amount of ${}^{\mathrm{14}}\mathrm{C}$ per gram of ice for ${}^{\mathrm{14}}\mathrm{CO}$, ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$, and ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ (Table 1) is calculated using the same method as in Petrenko et al. (2016) and is consistent with the Hippe and Lifton (2014) formulations for in situ ${}^{\mathrm{14}}\mathrm{C}$ concentrations. ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ and ${}^{\mathrm{14}}\mathrm{CO}$ measurements in our samples after all associated corrections, in addition to earlier Taylor Glacier results from Petrenko et al. (2016), are shown in Table 1 and Fig. 3.

2.3.2

Sublimation and processing of samples for ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements

CO${}_{\mathrm{2}}$ was liberated from ice samples using a newly developed ice sublimation device at the University of Rochester (Hmiel, 2020), roughly following the design of Schmitt et al. (2011). To briefly summarize the procedure, 1.5–2 kg ice samples were loaded into a vacuum glass vessel. The vessel was then evacuated, and the ice was sublimated at vacuum with six infrared emitters (Emitted Energy, USA) for 8–10 h. We did not sublimate 100 % of the samples because, as the ice sublimates away, impurities such as dust and organics start to accumulate on the surface. The aggregation of impurities on the sublimation front might enhance unwanted chemical reactions that produce extraneous carbon (Schmitt et al., 2011). Furthermore, towards the end of the extraction, the sublimation became less efficient as less surface area was available to absorb radiation. Approximately 1 kg of ice was sublimated in 8–10 h. However, the incomplete sublimation does not compromise the 100 % extraction efficiency, as all the gases trapped in the ice that is sublimated away are still released (Schmitt et al., 2011).

The liberated CO${}_{\mathrm{2}}$ was cryogenically trapped with liquid nitrogen, and the air was also cryogenically trapped with 5Å molecular sieve (Sigma-Aldrich, USA) under liquid nitrogen. After the sublimation was completed, the trapped CO${}_{\mathrm{2}}$ and air were expanded into separate volume-calibrated manometers in which pressure measurements were taken to calculate the [CO${}_{\mathrm{2}}$]. Finally, the isolated CO${}_{\mathrm{2}}$ was cryogenically transferred to and flame-sealed into a Pyrex glass ampule. The CO${}_{\mathrm{2}}$ was graphitized at ANSTO using the micro-furnaces (Yang and Smith, 2017), and the graphitized samples were measured for ${}^{\mathrm{14}}\mathrm{C}$ activity at the ANTARES AMS facility (Smith et al., 2010). One ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ sample (replicate for 30 m depth sample) was unfortunately lost during sublimation because the ice fractured under vacuum during the evacuation step.

A $\sim \mathrm{50}$–75 g ice subsample was taken from every ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ sample and shipped to OSU. The aliquots were measured for [CO${}_{\mathrm{2}}$], following Ahn et al. (2009), and [${\mathrm{CH}}_{\mathrm{4}}$] and TAC, following Mitchell et al. (2013; Table S8). Five field-produced bubble-free ice (BFI) samples and nine laboratory-produced BFI samples were also sublimated along with the glacial ice samples. During the sublimation of the BFI samples, a standard gas with known ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ activity and [CO${}_{\mathrm{2}}$] was introduced into the bottom of the glass sublimation vessel at 0.15 scc min${}^{-\mathrm{1}}$ flow rate for 8–10 h. The flow rate was set to mimic the rate of air liberation from glacial ice samples and the processing time also mimicked the amount of time needed to sublimate glacial ice samples. We used a standard gas with dead ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ activity for four laboratory-produced BFI samples and a standard gas with modern ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ activity for the other five laboratory-produced BFI samples. The CO${}_{\mathrm{2}}$ was cryogenically trapped downstream, processed, and measured for ${}^{\mathrm{14}}\mathrm{C}$ activity, following the same methods as the ice samples. In combination with the OSU [CO${}_{\mathrm{2}}$] and TAC measurements, the BFI samples were used to constrain the amount of extraneous carbon and ${}^{\mathrm{14}}\mathrm{C}$ introduced by sample transport, storage, and processing (see Sect. S1; Table S9). Finally, 11 commensurately sized ${}^{\mathrm{14}}\mathrm{C}$ standards and blanks (14–16 $\mathrm{\mu }$gC) with known ${}^{\mathrm{14}}\mathrm{C}$ activities (in 0–135 pMC range) were prepared, graphitized, and measured at ANSTO, concurrently with all the samples (Table S7), to characterize the effects from the addition of extraneous C during the graphitization process.

The detailed corrections for the ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ samples are discussed in the Supplement. We correct for the effects of extraneous C from graphitization and other ANSTO processing using a linear empirical correction from the commensurately sized ${}^{\mathrm{14}}\mathrm{C}$ standards (Sect. S1.1; Fig. S2b). The effects of extraneous carbon from ice sublimation/CO${}_{\mathrm{2}}$ extraction are calculated from the difference in measured ${}^{\mathrm{14}}\mathrm{C}$ activity of the laboratory-produced BFI samples relative to the measured ${}^{\mathrm{14}}\mathrm{C}$ activity of the standard gases with a mass balance approach (Sect. S1.2; Tables S8 and S9). Finally, the samples were corrected for the effects of post-coring in situ ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ production in ice using results from the field-produced BFI samples (Sect. S1.3; Table S10). The ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements in our samples after all associated corrections with their error-propagated uncertainties are shown in Table 1 and Fig. 3.

Figure 4

(a) ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}{/}^{\mathrm{14}}\mathrm{CO}$. (b) ${}^{\mathrm{14}}\mathrm{CO}$/total ${}^{\mathrm{14}}\mathrm{C}$ and ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$/total ${}^{\mathrm{14}}\mathrm{C}$ ratio. The solid lines represent the mean, and the dashed lines represent 2 standard deviations of the ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}{/}^{\mathrm{14}}\mathrm{CO}$ ratio, ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, and ${}^{\mathrm{14}}\mathrm{CO}$ fractions for samples deeper than 6.85 m where production by muons dominates. The ratios for rejected samples (Sects. 4.1, S2) are not shown. All error bars shown in this figure are 95 % CI.

[Figure omitted. See PDF]

An in-depth discussion about the analytical uncertainty in the ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements obtained from the sublimation method (which is important for the interpretation of the data because it is the largest source of uncertainty for total ${}^{\mathrm{14}}\mathrm{C}$) is provided in Sect. S1.4. In brief, we used the pooled standard deviation of replicate pairs ($\pm \mathrm{7.8}$ ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ molec. per gram of ice; $\mathrm{2}\mathit{\sigma }$) as the uncertainty for all ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements, except the 2.25 m sample pair (where we used the error-propagated uncertainties instead; Table 1). For the rest of the paper, we refer to the sum of measured ${}^{\mathrm{14}}\mathrm{CO}$, ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ as the total ${}^{\mathrm{14}}\mathrm{C}$. The ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}{/}^{\mathrm{14}}\mathrm{CO}$, ${}^{\mathrm{14}}\mathrm{CO}$/total and ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$/total ${}^{\mathrm{14}}\mathrm{C}$ ratio of the samples are shown in Fig. 4.

2.4 Sample integrity

Several samples were excluded from the data analysis; the detailed reasoning for rejecting these samples is discussed in Sects. S2 and S3. The surface samples (0.75 m depth) for all three ${}^{\mathrm{14}}\mathrm{C}$ species (${}^{\mathrm{14}}\mathrm{CO}$, ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$, and ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$) are rejected because of ambient air contamination from abundant fractures in the ice sample (due to thermal stresses in near-surface ice) and the likely chemical and/or biological alteration of ${\mathrm{CH}}_{\mathrm{4}}$, CO, and CO${}_{\mathrm{2}}$. The 19.5 and 30 m ${}^{\mathrm{14}}\mathrm{CO}$ and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ samples from the 2015/2016 season were rejected because of anomalous alterations in TAC, [CO], and [${\mathrm{CH}}_{\mathrm{4}}$]. The 30 m ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ sample was also rejected due to alteration in [CO${}_{\mathrm{2}}$].

3

Estimating the muogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates

3.1 Ice flow model to constrain sample exposure history

We used a 2D ice flow model from Buizert et al. (2012) to generate flow trajectories for the ice parcels corresponding to each sample depth (Fig. 5). The model first computes the 2D steady-state ice flow velocity field based on the observed surface velocities, ablation rates, and glacier thickness and then generates an ice parcel back-trajectory using a 2D linear interpolation of the ice flow velocity field (Kavanaugh et al., 2009a, b; Kavanaugh and Cuffey, 2009; Bliss et al., 2011). The largest source of uncertainty for the trajectories are the ablation rates (Buizert et al., 2012), which are based on measurements of 163 poles initially planted in 2002/2003 (Kavanaugh et al., 2009b; Bliss et al., 2011). All survey poles were measured a year later, providing 1-year average ablation rate estimates (Kavanaugh et al., 2009b) and again in the 2006/2007 season. Additionally, 17 poles were remeasured in the 2009/2010 and 2010/2011 seasons (Buizert et al., 2012). The 4-year average ablation rates were 4.7 cm yr${}^{-\mathrm{1}}$ higher than the 1-year average (Kavanaugh et al., 2009b); following Buizert et al. (2012), the ablation rate uncertainty for each pole was calculated by dividing 4.7 cm yr${}^{-\mathrm{1}}$ with $\sqrt{N}$, where $N$ is the length of the observation period in years ($N=\mathrm{1}$, 4, 7, or 8). Figure S3 shows the ablation rates along the glacier and their uncertainties inferred from the survey pole data.

The ice flow model used a bedrock profile from Kavanaugh et al. (2009a); however, the bedrock profile only extends to 72 km away from the glacier terminus, a point which we refer to as the glacier head (Fig. 5). This bedrock profile corresponds to 5–6 kyr of ice flow history and approximately one ${}^{\mathrm{14}}\mathrm{C}$ half-life. Beyond the constraints from the bedrock profile, we had to make assumptions about the depth of the long-term transport ($z_{\mathrm{deep}}$). Morse et al. (1998) provided a radar-based bedrock profile that includes the Taylor Glacier snow accumulation area (Baggenstos et al., 2018) north of the Taylor Dome – approximately 60 km upstream from the glacier head, where the Kavanaugh et al. (2009a) bedrock profile ends. Based on the bedrock profile from Morse et al. (1998), at the Taylor Glacier accumulation area, the depth of $\sim \mathrm{80}$ kyr ice (which corresponds to our 72 m sample) is $\sim \mathrm{575}$ m. We thus assumed that the depth of long-term transport ($z_{\mathrm{deep}}$) for the 72 m sample under the best-estimate ablation rate scenario (which we define as the reference sample) is 575 m. For other ice parcel trajectories ($i$), we scaled the depth of long-term transport ($z_{\mathrm{deep}}$) as follows:

5$$z_{\mathrm{deep}}\left(i\right)=\mathrm{575}-(z_{\mathrm{ref}}-z_{\mathrm{head}}(i\left)\right),$$ where $z_{\mathrm{ref}}$ represents the depth of the 72 m reference sample in the model at the glacier head under the best-estimate flowline ($z_{\mathrm{ref}}$ is 699 m), and $z_{\mathrm{head}}$ represents the depth of the ice parcel of interest at the glacier head. We assumed that the difference in depth between the reference sample and the sample of interest ($i$) at the glacier head and during long-term transport within the glacier is the same.

Figure 5

Example of the ice parcel back-trajectory and associated uncertainties. For the Monte Carlo estimate of uncertainties (Sect. 5.4.2), for each given sample depth (72 m in this figure), 10 000 back-trajectories are generated. Each back-trajectory corresponds to a different ablation rate scenario (the ablation rates are perturbed within their experimental measurement uncertainties to generate the scenarios). The shaded region represents the 68 % CI uncertainty envelope of the flow trajectory.

[Figure omitted. See PDF]

Figure 6

(a) Comparison between total ${}^{\mathrm{14}}\mathrm{C}$ measurements with modeled best-estimate ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ parameters from this study and Heisinger et al. (2002a, b). (b) Comparison between total ${}^{\mathrm{14}}\mathrm{C}$ measurements with modeled best-estimate ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ parameters from this study and modeled total ${}^{\mathrm{14}}\mathrm{C}$ with best-fit $f_{\mathrm{tot}}$ when ${\mathit{\sigma }}_{\mathrm{0}}$ is forced to be zero. (c) Comparison between ${}^{\mathrm{14}}\mathrm{CO}$ measurements with modeled best-estimate $f_{\mathrm{neg}}$ and $f_{\mathrm{fast}}$ parameters from this study and modeled ${}^{\mathrm{14}}\mathrm{CO}$ with best-fit $f_{\mathrm{neg}}$ when $f_{\mathrm{fast}}$ is forced to be zero. (d) Comparison between total ${}^{\mathrm{14}}\mathrm{C}$ measurements with modeled best-estimate ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ parameters from this study and modeled total ${}^{\mathrm{14}}\mathrm{C}$ from the sensitivity analyses when we assume 25 % contribution from organics and high ablation rate scenario (Fig. S9). The thin colored lines represent the 95 % CI envelope of the model results (corresponding to the contour plot of Fig. 7a for Fig. 6a and b and the contour plot of Fig. 7b for Fig. 6c). The error bars shown in the data are 95 % CI. In Fig. 6d, the solid blue line represents the ${}^{\mathrm{14}}\mathrm{C}$ profile from the modeled best-estimate $f_{\mathrm{tot}}$ and ${\mathit{\sigma }}_{\mathrm{0}}$ under the best-estimate ablation rate. The solid red line represents the ${}^{\mathrm{14}}\mathrm{C}$ profile from the sensitivity analysis when $f_{\mathrm{tot}}$ and ${\mathit{\sigma }}_{\mathrm{0}}$ are tuned to fit the total ${}^{\mathrm{14}}\mathrm{C}$ data that are scaled by 25 % to account for contribution from organics (red dots) under high ablation rate scenario (Fig. S9). The dashed red line represents the ${}^{\mathrm{14}}\mathrm{C}$ profile from the sensitivity analysis when $f_{\mathrm{tot}}$ is kept constant at 0.021 (which is the minimum value provided by Heisinger et al., 2002b) and ${\mathit{\sigma }}_{\mathrm{0}}=\mathrm{0}$ (which provides the best-fit against the data).

[Figure omitted. See PDF]

3.2

${}^{\mathrm{14}}\mathrm{C}$ production in sample ice parcel

We used the model for in situ cosmogenic nuclide production by muons from Balco et al. (2008) and Balco (2017; model 1A in Balco, 2017), with all relevant parameters adjusted for ice (Fig. 1). This model in turn uses the Heisinger et al. (2002a, b) parameterizations described above and the additional altitude scaling of the muon fluxes described in Balco et al. (2008). We then used a forward model that numerically integrates the total ${}^{\mathrm{14}}\mathrm{C}$ in the ice sample along its flow path. For the initial condition, we assumed that, at the depth of long-term transport ($z_{\mathrm{deep}}$), the ${}^{\mathrm{14}}\mathrm{C}$ concentration in the ice parcel is at the steady state, as follows: 6$${\displaystyle \frac{\mathrm{d}C}{\mathrm{d}t}}\left(\mathrm{at}{z}_{\mathrm{deep}}\right)=\mathrm{0}=P_{\mathrm{neg}}\left(z_{\mathrm{deep}}\right)+P_{\mathrm{fast}}\left(z_{\mathrm{deep}}\right)-C_{\mathrm{0}}\mathit{\lambda }.$$ The steady-state assumption means that, at $z_{\mathrm{deep}}$, the rate of radioactive decay ($C_{\mathrm{0}}\mathit{\lambda }$) is balanced by production from negative muon capture ($P_{\mathrm{neg}}$) and fast muon reaction ($P_{\mathrm{fast}}$). We use a ${}^{\mathrm{14}}\mathrm{C}$ decay constant $\mathit{\lambda }=\mathrm{1.216}\times {\mathrm{10}}^{-\mathrm{4}}$ yr${}^{-\mathrm{1}}$ for all our calculations (including the conversion of ${}^{\mathrm{14}}\mathrm{C}$ units in the previous section), which corresponds to the recommended ${}^{\mathrm{14}}\mathrm{C}$ half-life of 5700 years (Kutschera, 2019). For each ice parcel, we calculated the steady state, the initial ${}^{\mathrm{14}}\mathrm{C}$ concentration ($C_{\mathrm{0}}$) from Eq. (6), and then used the following differential equation, 7$${\displaystyle \frac{\mathrm{d}C}{\mathrm{d}t}}=P_{\mathrm{neg}}\left(z\left(t\right)\right)+P_{\mathrm{fast}}\left(z\left(t\right)\right)-C\mathit{\lambda },$$ to numerically integrate the ${}^{\mathrm{14}}\mathrm{C}$ concentration of the ice parcel along the flow trajectory. To avoid interference from spallogenic ${}^{\mathrm{14}}\mathrm{C}$, we only considered samples with mean depths greater than 6.85 m depth.

We sampled the parameter space in a grid search approach to obtain the best-estimate values for muogenic ${}^{\mathrm{14}}\mathrm{C}$ production parameters ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$. Using the best-estimate flow trajectory, we calculated the expected ${}^{\mathrm{14}}\mathrm{C}$ in the samples corresponding to all combinations of ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$, with each of the parameters ranging between 0 % and 100 % of the values from Heisinger et al. (2002a, b). To save computational time, we first conducted the grid search at a coarse resolution of 10 % increments (Fig. S4a). The goodness of the fit ($\mathit{\chi }\mathrm{2}$) for each simulation was calculated as follows: 8$${\mathit{\chi }}^{\mathrm{2}}=\sum {\displaystyle \frac{(C_{\mathrm{obs}\left(z\right)}-C_{\exp }(z){)}^{\mathrm{2}}}{C_{\exp }\left(z\right)}},$$ where $C_{\mathrm{obs}}\left(z\right)$ is the measured total ${}^{\mathrm{14}}\mathrm{C}$, and $C_{\exp }\left(z\right)$ is the total ${}^{\mathrm{14}}\mathrm{C}$ (${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ $+$ ${}^{\mathrm{14}}\mathrm{CO}$ $+$ ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$; Fig. 3d) calculated by the forward model at sample depth $z$. To find more precise best-estimate ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$, we conducted the grid search again at a higher resolution in 0.2 % increments from the Heisinger et al. (2002a, b) values near the $\mathit{\chi }\mathrm{2}$ minimum, between 0 and 0.0352 mb for ${\mathit{\sigma }}_{\mathrm{0}}$ and 0 to 0.01 for $f_{\mathrm{tot}}$ (Fig. S4b).

To estimate the uncertainties in ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$, we used a Monte Carlo sampling of model parameters. We assumed that the ablation rate uncertainties (Fig. S3) represent $\mathrm{2}\mathit{\sigma }$ normally distributed uncertainties. We then perturbed the ablation rates within their uncertainties and generated a pool of 10 000 possible flow trajectories for each sample depth. However, in 69 out of 10 000 flow scenarios, the ice parcel back-trajectories hit the bedrock and became unphysical afterwards. These unphysical trajectories were removed from the pool of possible ice flow trajectories. Next, we started with the best-estimate ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ and assumed a normally distributed and large 200 % ($\mathrm{1}\mathit{\sigma }$) error for each parameter (Fig. S5a) as the prior distribution for the Monte Carlo method. We removed ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ values that are below zero from the prior distribution because they are unphysical and conducted 100 000 Monte Carlo simulations using the forward ${}^{\mathrm{14}}\mathrm{C}$ production model. For each Monte Carlo simulation, we randomly picked one of the previously generated possible ice flow trajectories and a random pair of ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ from the generated prior distributions (Fig. S5a). We then calculated the expected ${}^{\mathrm{14}}\mathrm{C}$ concentrations for each sample depth using the forward model and compared the model–data fit. We accept all pairs of ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ values that produce the model-calculated total ${}^{\mathrm{14}}\mathrm{C}$ within the 95 % CI (7.8 ${}^{\mathrm{14}}\mathrm{C}$ atoms per gram of ice) and 67 % CI (3.9 ${}^{\mathrm{14}}\mathrm{C}$ atoms per gram of ice) analytical uncertainty of the best-fit, model-calculated total ${}^{\mathrm{14}}\mathrm{C}$ (black line in Fig. 6). The ranges of accepted ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ pairs are shown in Fig. 7a as contours. The discussion about the selection of acceptance criteria for estimating ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ uncertainties is provided in Sect. S1.4.

Figure 7

(a) The 68 % and 95 % CI contours of accepted ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ values for total ${}^{\mathrm{14}}\mathrm{C}$. (b) The 68 % and 95 % CI contours of accepted $f_{\mathrm{neg}}$ and $f_{\mathrm{fast}}$ values for ${}^{\mathrm{14}}\mathrm{CO}$ (see Sect. 5.4.3). For comparison, the ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ values from Heisinger et al. (2002a, b) are shown as a blue star. The best-fit values for ${\mathit{\sigma }}_{\mathrm{0}}$, $f_{\mathrm{tot}}$, $f_{\mathrm{neg}}$, and $f_{\mathrm{fast}}$ are shown as a red star in both figures. The best-fit values for ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ are when $\mathit{\alpha }$ (energy scaling in Eq. 3).

[Figure omitted. See PDF]

3.3

${}^{\mathrm{14}}\mathrm{CO}$ production model in sample ice parcel

The in situ cosmogenic ${}^{\mathrm{14}}\mathrm{CO}$ production rates in ice are of specific interest, as discussed in Sect. 1.1. To characterize the ${}^{\mathrm{14}}\mathrm{CO}$ production rates, we introduced additional scaling factors $f_{\mathrm{neg}}$ and $f_{\mathrm{fast}}$ for negative muon and fast muon mechanisms, respectively, as tuned model parameters. The differential equation of Eq. (7) is modified into the following: 9$${\displaystyle \frac{\mathrm{d}{(}^{\mathrm{14}}\mathrm{CO})}{\mathrm{d}t}}=f_{\mathrm{neg}}{P}_{\mathrm{neg}}\left(z\left(t\right)\right)+f_{\mathrm{fast}}{P}_{\mathrm{fast}}\left(z\left(t\right)\right)-{(}^{\mathrm{14}}\mathrm{CO})\mathit{\lambda }.$$ We note that $P_{\mathrm{neg}}$ and $P_{\mathrm{fast}}$ in Eq. (9) are the total ${}^{\mathrm{14}}\mathrm{C}$ production rates calculated from the Balco et al. (2008) model. The scaling factors $f_{\mathrm{neg}}$ and $f_{\mathrm{fast}}$ each encompass two terms, i.e., one that adjusts the total ${}^{\mathrm{14}}\mathrm{C}$ production rates and another that accounts for the ${}^{\mathrm{14}}\mathrm{CO}$ fraction of total ${}^{\mathrm{14}}\mathrm{C}$. The determination of best-estimate $f_{\mathrm{neg}}$ and $f_{\mathrm{fast}}$ and their uncertainties was similar to the approach for ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ described above. The $\mathit{\chi }\mathrm{2}$ grid search was conducted with all combinations of $f_{\mathrm{neg}}$ and $f_{\mathrm{fast}}$ values ranging from 0 to 0.2 at 0.001 resolution (Fig. S4c). Similar to the total ${}^{\mathrm{14}}\mathrm{C}$ data, we used the average analytical uncertainty in the ${}^{\mathrm{14}}\mathrm{CO}$ sample set as the acceptance criteria for the Monte Carlo simulations to estimate the uncertainties of $f_{\mathrm{neg}}$ and $f_{\mathrm{fast}}$. We accepted all sets of $f_{\mathrm{neg}}$ and $f_{\mathrm{fast}}$ from the 100 000 Monte Carlo simulations that yielded model-predicted ${}^{\mathrm{14}}\mathrm{CO}$ within 1.2 ${}^{\mathrm{14}}\mathrm{CO}$ molec. per gram of ice (95 % CI uncertainty) and 0.6 ${}^{\mathrm{14}}\mathrm{CO}$ molec. per gram of ice (68 % CI uncertainty) from the best-fit model (Fig. 8). Figure 7b shows (as contours) the accepted sets of $f_{\mathrm{neg}}$ and $f_{\mathrm{fast}}$ values.

Figure 8

(a) Comparison between measured total ${}^{\mathrm{14}}\mathrm{C}$ from Scharffenbergbotnen, expected total ${}^{\mathrm{14}}\mathrm{C}$ using the production rates inferred in this study, and expected total ${}^{\mathrm{14}}\mathrm{C}$ using the Heisinger et al. (2002a, b) production rates. (b) Comparison between measured ${}^{\mathrm{14}}\mathrm{CO}$ from Scharffenbergbotnen and expected ${}^{\mathrm{14}}\mathrm{CO}$, using the production rates inferred from Taylor Glacier. The colored lines in both figures represent the 95 % CI envelope of the model results. At the depths plotted in this figure (deeper than 5 m), the production from neutron-induced spallation is negligible.

[Figure omitted. See PDF]

3.4 Comparison with Scharffenbergbotnen ablation site

Van der Kemp et al. (2002) measured ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ and ${}^{\mathrm{14}}\mathrm{CO}$ in ice from the Scharffenbergbotnen ice ablation site, Antarctica. Using a 1D ablation model, we examined how the estimates of muogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates from Taylor Glacier compare to the Scharffenbergbotnen data. We assumed that the measured ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ plus ${}^{\mathrm{14}}\mathrm{CO}$ from Scharffenbergbotnen are comparable to our measurements of total ${}^{\mathrm{14}}\mathrm{C}$ in Taylor Glacier (since our data show that less than 0.3 % of the total ${}^{\mathrm{14}}\mathrm{C}$ from muon production forms ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$; Sect. 4.1). We then used the ${}^{\mathrm{14}}\mathrm{C}$ concentration from the deepest Scharffenbergbotnen sample (45 m) as the initial condition. In the 1D ablation model, the Scharffenbergbotnen ice parcel moves upward at a rate ($\mathrm{d}z/\mathrm{d}t$) equal to the ablation rate from stake measurements (Eq. 10; $a=\mathrm{16}\pm \mathrm{4}$ cm yr${}^{-\mathrm{1}}$). $$\begin{array}{}\text{10}& {\displaystyle \frac{\mathrm{d}z}{\mathrm{d}t}}=-a\text{11}& {\displaystyle \frac{\mathrm{d}C}{\mathrm{d}t}}=P_n^{\prime }\left(z\left(t\right)\right)+{P^{\prime }}_{\mathrm{neg}}\left(z\left(t\right)\right)+{P^{\prime }}_{\mathrm{fast}}\left(z\left(t\right)\right)-C\mathit{\lambda }.\end{array}$$

The expected ${}^{\mathrm{14}}\mathrm{C}$ concentration in the ice is given by the differential equation (Eq. 11), where $P_n^{\prime }$ is the ${}^{\mathrm{14}}\mathrm{C}$ spallogenic production rate from Young et al. (2014; $\mathrm{12.0}\pm \mathrm{0.9}$ atoms per gram of quartz per year at the surface), first scaled to the SLHL production rate in ice ($\mathrm{20.0}\pm \mathrm{1.5}$ atoms per gram of ice per year at the surface) accounting for the number of ${}^{\mathrm{16}}\mathrm{O}$ atoms per gram of ice vs. quartz (variable $N$; Eq. 3; Petrenko et al., 2016) and then to the production rate at the Scharffenbergbotnen site (1173 m above sea level or m a.s.l.) using the Lifton et al. (2014) LSDn nuclide-specific model ($P_n^{\prime }=\mathrm{71.2}\pm \mathrm{3.6}$ atoms per gram of ice per year at the surface). $P_{\mathrm{neg}}^{\prime }$ and $P_{\mathrm{fast}}^{\prime }$ are the muogenic production rates inferred from the Taylor Glacier data scaled to the elevation of Scharffenbergbotnen (1173 m a.s.l.), using the altitude scaling factors from Balco (2017). We also repeated this calculation for ${}^{\mathrm{14}}\mathrm{CO}$ only to compare the muogenic ${}^{\mathrm{14}}\mathrm{CO}$ production rates with the ${}^{\mathrm{14}}\mathrm{CO}$ data from Scharffenbergbotnen.

4 Results and discussion

4.1

Measured ${}^{\mathrm{14}}\mathrm{C}$ values and partitioning of ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, ${}^{\mathrm{14}}\mathrm{CO}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$

Table 1 and Fig. 3a–c show the depth profiles of ${}^{\mathrm{14}}\mathrm{CO}$, ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$, and ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ after all corrections. For the ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements, comparison with prior results that used a wet extraction approach (Fig. S6) confirms the caveats discussed by Petrenko et al. (2016) in that their ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements were uncertain and represent the upper bound. The ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}{/}^{\mathrm{14}}\mathrm{CO}$ ratios from the new samples ($\mathrm{0.0074}\pm \mathrm{0.0004}$; 95 % CI; $n=\mathrm{4}$; from all samples below 19.5 m) appear to be constant within uncertainties (Fig. 4a), which is in agreement with earlier results ($\mathrm{0.0076}\pm \mathrm{0.0004}$; 95 % CI; $n=\mathrm{4}$) from Petrenko et al. (2016). This confirms that the two muon reactions produce ${}^{\mathrm{14}}\mathrm{C}$ in a constant ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}{/}^{\mathrm{14}}\mathrm{CO}$ ratio. The ${}^{\mathrm{14}}\mathrm{CO}$ and ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ fractions of total ${}^{\mathrm{14}}\mathrm{C}$ are also relatively constant at depth (Fig. 4b), suggesting that the two muon reactions produce all three ${}^{\mathrm{14}}\mathrm{C}$ species in constant ratios.

For samples deeper than 6.85 m, on average 33.7 % ($\pm \mathrm{11.4}\mathit{\%}$; 95 % CI) of the produced cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ becomes ${}^{\mathrm{14}}\mathrm{CO}$ and 66.1 % ($\pm \mathrm{11.5}\mathit{\%}$; 95 % CI) of the produced cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ becomes ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ (Fig. b). The uncertainties of ${}^{\mathrm{14}}\mathrm{CO}$ and ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ fractions on the deepest samples (72 m depth) are relatively large because of the small ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ signal (11.8 to 13.6 ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ molec. per gram of ice) relative to the uncertainty in our measurements ($\pm \mathrm{7.8}$ ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ molec. per gram of ice; 95 % CI). The ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ fraction in samples that are deeper than 6.85 m ($\mathrm{0.66}\pm \mathrm{0.12}$; 95 % CI) is also in agreement with prior reported ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ fraction of 0.69 from the Scharffenbergbotnen ablation site (Van der Kemp et al., 2002). Finally, the shallow samples ($<\mathrm{6}$ m ice equivalent) show higher ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}{/}^{\mathrm{14}}\mathrm{CO}$ ratios (Fig. 4a) and ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ $/$ total ${}^{\mathrm{14}}\mathrm{C}$ ratios. This may indicate that spallation produces higher amounts of ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ and ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ relative to ${}^{\mathrm{14}}\mathrm{CO}$ (Petrenko et al., 2016) or that CO (and ${}^{\mathrm{14}}\mathrm{CO}$) is not well-preserved in near-surface ice of the Taylor Glacier.

At depths where production from muons dominates ($>\mathrm{6}$ m), less than 0.3 % of the produced cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ in ice forms ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ (Table 1; Fig. 4). Although the ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ measurement from the 10 m depth sample is not available (Petrenko et al., 2016), we still include the 10 m data point in the total ${}^{\mathrm{14}}\mathrm{C}$ dataset used to infer ${\mathit{\sigma }}_{\mathrm{0}}$ and $f_{\mathrm{tot}}$ values and their uncertainties. The contribution from ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ (which would have been on the order of $\sim \mathrm{0.2}$ ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ molec. per gram of ice; Fig. 3b) is insignificant compared to the uncertainty in total ${}^{\mathrm{14}}\mathrm{C}$. We account for the lack of ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ measurement at this depth by scaling the total ${}^{\mathrm{14}}\mathrm{C}$ of the 10 m sample by a factor of $\mathrm{1.003}\pm \mathrm{0.003}$ (95 % CI; Table 1).

Table 2

Probability of the ${}^{\mathrm{14}}\mathrm{C}$ production from stopped negative muons ($f_{\mathrm{tot}}$), reference nuclear reaction cross section for production via fast muon interactions (${\mathit{\sigma }}_{\mathrm{0}}$), and total ${}^{\mathrm{14}}\mathrm{C}$ production rates in ice at the surface from the two muon reactions rescaled to SLHL (sea level high latitude), using the Lifton et al. (2014) LSDn scaling. All errors shown represent the 95 % confidence interval.

	Overall probability of	Reference nuclear	SLHL total ${}^{\mathrm{14}}\mathrm{C}$	SLHL total ${}^{\mathrm{14}}\mathrm{C}$
	negative muon capture	reaction cross section	production rate in ice	production rate in ice
	reaction ($f_{\mathrm{tot}}$)	(${\mathit{\sigma }}_{\mathrm{0}}$; mb)	by negative muons	by fast muons
			(atoms per gram of ice per year)	(atoms per gram of ice per year)
This study	0.0044 ($+\mathrm{0.0026}/-\mathrm{0.0026}$)	0.0024 ($+\mathrm{0.0017}/-\mathrm{0.0018}$)	0.79 ($+\mathrm{0.47}/-\mathrm{0.46}$)	0.21 ($+\mathrm{0.16}/-\mathrm{0.15}$)
Heisinger et al. (2002a, b)	$\mathrm{0.025}\pm \mathrm{0.004}$	0.0088 ($+\mathrm{0.0098}/-\mathrm{0.0088}$)	$\mathrm{4.76}\pm \mathrm{0.76}$	0.74 ($+\mathrm{0.83}/-\mathrm{0.74}$)
Lupker et al. (2015)	0.024 ($+\mathrm{0.006}/-\mathrm{0.016}$)${}^{*}$	0 ($+\mathrm{0.0118}/-\mathrm{0}$)	4.70 ($+\mathrm{1.22}/-\mathrm{3.04}$)	0 ($+\mathrm{1.52}/-\mathrm{0}$)

${}^{*}$ Adjusted to ice, assuming that the chemical compound factor ($f_{\mathrm{c}}$) of ice is 1.0, and $f_{\mathrm{c}}$ for quartz is 0.704 (Heisinger et al., 2002b).

4.2

Inferred muogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates in ice and comparison with production rates in quartz

Assuming that the majority of in situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ in ice forms ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, ${}^{\mathrm{14}}\mathrm{CO}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$, the muogenic ${}^{\mathrm{14}}\mathrm{C}$ production parameters from Heisinger et al. (2002a, b; $f_{\mathrm{tot}}$ for negative muon capture and ${\mathit{\sigma }}_{\mathrm{0}}$ for fast muon reaction) are well outside the confidence intervals of our measurements (Table 2; Figs. 6a and 7a). Using the larger uncertainty for ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements obtained from step-by-step error propagation (Sect. S1.4; Fig. S7) does not change this conclusion. We calculated factors of 5.7 (3.6–13.9; 95 % CI) and 3.7 (2.0–11.9; 95 % CI) lower probability of the negative muon capture reaction ($f_{\mathrm{tot}}$) and the reference cross section for the fast muon mechanism (${\mathit{\sigma }}_{\mathrm{0}}$) compared to the values given by Heisinger et al. (2002a, b).

One possible explanation for the disagreement between our results and those of Heisinger et al. (2002a, b) is that our ${}^{\mathrm{14}}\mathrm{C}$ measurements (mostly either ${}^{\mathrm{14}}\mathrm{CO}$ or ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, as ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ only constitutes $<\mathrm{0.3}\mathit{\%}$ of total ${}^{\mathrm{14}}\mathrm{C}$) might be incorrect. However, in the following we thoroughly explore this possibility and argue that it is very unlikely. Our ${}^{\mathrm{14}}\mathrm{CO}$ measurements used a well-established analytical technique (e.g., Dyonisius et al., 2020; Hmiel et al., 2020; Petrenko et al., 2013, 2017, 2021). With regards to ${}^{\mathrm{14}}\mathrm{CO}$ measurements in air, this analytical technique (Petrenko et al., 2021) yields comparable results to independent, atmospheric ${}^{\mathrm{14}}\mathrm{CO}$ measurements from other research groups (e.g., Manning et al., 2005; Mak and Southon, 1998). We also have no reason to believe that there is a systematic loss of ${}^{\mathrm{14}}\mathrm{CO}$ during the ice melting process. The ice melting (wet extraction) ensures that all ${}^{\mathrm{14}}\mathrm{C}$ and CO are liberated. CO is not very soluble in water (the dissolved CO fraction at equilibrium in our system is on the order of 1 %), and we used the measured $\mathit{\delta }\mathrm{Xe}/\mathrm{Kr}$ (Table S4) to correct for the solubility effects for both [CO] and [${\mathrm{CH}}_{\mathrm{4}}$]. The on-site field extraction within hours of the sample retrieval ensures that there is minimal post-coring gas loss. Finally, ice core and firn air ${}^{\mathrm{14}}\mathrm{CO}$ measurements at Summit, Greenland, are consistent within uncertainties with muogenic ${}^{\mathrm{14}}\mathrm{CO}$ production rates inferred from Taylor Glacier (Hmiel et al., 2020; Hmiel, 2020).

The sublimation technique used for our ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements also ensures the 100 % extraction efficiency of gases (Schmitt et al., 2011). We tested the system performance by sublimating BFI (bubble-free ice) samples while adding standard gases with known ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ activities. As mentioned in Sect. 2.3.2, we used two standard gases with known ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ activities, i.e., one with modern ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ activity and the other with dead ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ activity, when sublimating the BFI samples. There is neither a significant alteration in the CO${}_{\mathrm{2}}$ mole fraction and ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ activity of both standard gases (Table S8) nor in the CO${}_{\mathrm{2}}$ mole fraction of the ice samples (Table S9), which suggests that the processes of sublimating ice and flowing gas through the system components do not result in loss of ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$. Finally, we can also rule out the possibility of post-coring ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ loss. In a separate measurement campaign (Hmiel, 2020), we brought the sublimation system to Summit, Greenland, and sublimated the ice samples on site (within days from the time when the ice core was drilled). We compared the ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ from the on-site field sublimation with the ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ from depth-adjacent replicates sublimated at the University of Rochester laboratory and found that they are indistinguishable within uncertainty.

We argue that another strong indication that our measurements are robust is the good agreement with independent results from Van der Kemp et al. (2002). Van der Kemp et al. (2002) measured ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ and ${}^{\mathrm{14}}\mathrm{CO}$ in ice from the Scharffenbergbotnen ice ablation site using a dry extraction technique. The total measured ${}^{\mathrm{14}}\mathrm{C}$ values were significantly lower than the expected values, based on the stake-measured ablation rates and muogenic production rates based on laboratory irradiations of quartz targets (Heisinger et al., 2000a, b). Van der Kemp et al. (2002) initially hypothesized that the low extraction efficiency of dry mechanical extraction (which can result in an incomplete release of the in situ produced ${}^{\mathrm{14}}\mathrm{C}$ from the ice grains) might be responsible for this discrepancy. However, we used a sublimation method for our ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements and a melt extraction method for our ${}^{\mathrm{14}}\mathrm{CO}$ measurements; both methods guarantee that all in situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ in the ice lattice is released. Figure 8 shows that the Scharffenbergbotnen data are consistent with the expected total ${}^{\mathrm{14}}\mathrm{C}$ and ${}^{\mathrm{14}}\mathrm{CO}$ from the Taylor-Glacier-derived production rates.

The good agreement in the ratio of ${}^{\mathrm{14}}\mathrm{C}$ compounds (${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ fraction is $\mathrm{0.66}\pm \mathrm{0.12}$ in this study and 0.69 in Van der Kemp et al., 2002) suggests that our extraction methods and analytical techniques were not systematically losing ${}^{\mathrm{14}}\mathrm{CO}$ or ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ (which would then bias the ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ and ${}^{\mathrm{14}}\mathrm{CO}$ fraction). It is theoretically possible that both our measurements and Van der Kemp et al. (2002) are wrong. However, to produce the same ratio of ${}^{\mathrm{14}}\mathrm{C}$ compounds, it would require all three analytical systems from these studies to be systematically wrong in the same direction and by the same magnitude, which is highly unlikely. The good agreement between the Taylor Glacier and Scharffenbergbotnen data suggests that dry mechanical extraction used by Van der Kemp et al. (2002) is a valid technique for extracting ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ and ${}^{\mathrm{14}}\mathrm{CO}$ from bubbly, non-clathrated ice cores. One possible explanation is that, after production, in situ ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ and ${}^{\mathrm{14}}\mathrm{CO}$ quickly migrates from the ice matrix to the air bubbles. This result is consistent with previous observations that the retention of in situ cosmogenic ${}^{\mathrm{14}}\mathrm{C}$ in firn grains is very low (Petrenko et al., 2013; Van der Kemp et al., 2000; Wilson and Donahue, 1990).

One way to reconcile our measurements and those of Van der Kemp et al. (2002) with the Heisinger et al. (2002a, b) production rates would be to have much higher (factor of 3 or greater) long-term (hundreds to thousands of years) ablation rates at Taylor Glacier and Scharffenbergbotnen compared to recent ablation stakes measurements. The ablation rate over the last 100 years is especially important with regards to in situ ${}^{\mathrm{14}}\mathrm{C}$ production rate from negative muon capture (Fig. S9). The ablation rate at blue ice areas is controlled by climate via a combination of temperature, insolation, and wind (mainly katabatic; e.g., Bintanja, 1999). To obtain much higher long-term ablation rates at both Taylor Glacier and Scharffenbergbotnen (which are on opposite sides of Antarctica), we would need either the temperatures to have dropped sharply in the last couple of decades, the winds to have slowed dramatically, or for insolation to have sharply decreased. However, the two deep ice cores nearest to Taylor Glacier, namely Taylor Dome (Steig et al., 2000) and RICE (Roosevelt Island Climate Evaluation; Bertler et al., 2018) have not shown large climate changes in the region over the last 1000 years. The European Project for Ice Coring in Antarctica (EPICA) Dronning Maud Land (EDML) ice core drilled nearby Scharffenbergbotnen blue ice area also has not shown large climate variability in the region over the Holocene period (EPICA Community Members, 2010). The glaciological survey of Taylor Glacier also indicated that the glacier is approximately at steady state, given the stake-measured ablation rates (Kavanaugh et al., 2009a, b). The 30-year record of weather observations over the McMurdo Dry Valleys area show that the recent climate in this region has been stable (Obryk et al., 2020). Finally, a 14-year-long observation study (Sinisalo et al., 2003) over the Scharffenbergbotnen blue ice area has also shown no significant recent change in the ablation rates. We thus argue that a large decrease in ablation rates in recent years, as compared to the long-term average (over the last hundreds to thousands of years), does not seem to be a realistic explanation.

Table 3

${}^{\mathrm{14}}\mathrm{CO}$-specific surface production rates in ice from the two muon mechanisms normalized to the SLHL (sea level high latitude) site using the Lifton et al. (2014) LSDn scaling. All errors shown represent the 95 % confidence interval.

	SLHL ${}^{\mathrm{14}}\mathrm{CO}$ production rate in ice by	SLHL ${}^{\mathrm{14}}\mathrm{CO}$ production rate in ice by
	negative muons (molec. per gram of ice per year)	fast muons (molec. per gram of ice per year)
This study	0.310 ($+\mathrm{0.075}/-\mathrm{0.063}$)	0.063 ($+\mathrm{0.022}/-\mathrm{0.018}$)
Petrenko et al. (2016)	0.24 ($+\mathrm{0.14}/-\mathrm{0.14}$)	0.053 ($+\mathrm{0.028}/-\mathrm{0.028}$)

For direct comparison with other studies, we used the scaling factors from the Lifton et al. (2014) LSDn model to calculate the corresponding sea level high latitude (SLHL) total ${}^{\mathrm{14}}\mathrm{C}$ and ${}^{\mathrm{14}}\mathrm{CO}$-specific production rates in ice (Tables 2 and 3). Our estimates of the ${}^{\mathrm{14}}\mathrm{CO}$-specific production rates agree with those of Petrenko et al. (2016), within errors (Table 3). Compared with the results from Petrenko et al. (2016), we also calculated a slightly smaller uncertainty in the ${}^{\mathrm{14}}\mathrm{CO}$-specific production rate by negative muon capture (Table 3). We converted the Lupker et al. (2015) estimates of $f_{\mathrm{tot}}$ in quartz into $f_{\mathrm{tot}}$ for ice (Table 2), using the chemical compound factors ($f_{\mathrm{C}}$) for quartz and ice from Heisinger et al. (2002b). With regards to the negative muon capture, the Lupker et al. (2015) estimate of $f_{\mathrm{tot}}$ is in close agreement with Heisinger et al. (2002b; Table 2). This result is supported by Balco (2017), who refitted the Leymon High data and obtained similar $f_{\mathrm{tot}}$ estimates. However, the high $f_{\mathrm{tot}}$ in Lupker et al. (2015), as compared to our result, was offset by their best ${\mathit{\sigma }}_{\mathrm{0}}$ estimate of zero (lower than our result). For a direct comparison with results from Lupker et al. (2015), we fit our data, while forcing ${\mathit{\sigma }}_{\mathrm{0}}$ (and hence ${}^{\mathrm{14}}\mathrm{C}$ production from fast muons) to be zero (Fig. 6b) and cannot find a scenario with reasonable model–data agreement.

Because of the relatively large uncertainty in the ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements, the total ${}^{\mathrm{14}}\mathrm{C}$ data still allow ${\mathit{\sigma }}_{\mathrm{0}}$ to be close to zero, given a sufficiently large $f_{\mathrm{tot}}$ (Fig. 7a). However, our ${}^{\mathrm{14}}\mathrm{CO}$ data (which have much lower relative uncertainties and use a more established measurement technique) unambiguously show that ${\mathit{\sigma }}_{\mathrm{0}}$ and ${}^{\mathrm{14}}\mathrm{C}$ production from fast muons cannot be zero (Figs. 8, 7b). As discussed in Lupker et al. (2015) and Balco (2017), the ${}^{\mathrm{14}}\mathrm{C}$ data from the 15.5 m Leymon High quartzite core might not cover the depth range in which production from fast muons dominates. In contrast, when integrated over the whole flow history, production from fast muons represents the dominant source of ${}^{\mathrm{14}}\mathrm{C}$ in our deeper samples. The very high-end estimate of our reference nuclear reaction cross section ${\mathit{\sigma }}_{\mathrm{0}}$ (for ${}^{\mathrm{14}}\mathrm{C}$ production from fast muons) is still within the large uncertainty of ${\mathit{\sigma }}_{\mathrm{0}}$ from both Heisinger et al. (2002a) and Lupker et al. (2015; Table 2). However, our estimated total probability of negative muon capture ($f_{\mathrm{tot}}$; and hence the ${}^{\mathrm{14}}\mathrm{C}$ production rate from negative muon) is well outside the confidence intervals of $f_{\mathrm{tot}}$ reported by both Heisinger et al. (2002b) and Lupker et al. (2015; Table 2).

One caveat to our estimated in situ muogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates in ice (and that of Van der Kemp et al., 2002) is that the total ${}^{\mathrm{14}}\mathrm{C}$ from the gas species we measured (${}^{\mathrm{14}}\mathrm{CO}$, ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$) might not account for all the muogenic in situ ${}^{\mathrm{14}}\mathrm{C}$. Although ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ and ${}^{\mathrm{14}}\mathrm{CO}$ likely constitute the large majority (Lal et al., 1997, 2000), a small amount of in situ ${}^{\mathrm{14}}\mathrm{C}$ can also form ${}^{\mathrm{14}}\mathrm{C}$-bearing organic materials. Measurements of ${}^{\mathrm{14}}\mathrm{C}$ in organic carbon from alpine ice for the purpose of radiocarbon dating have shown elevated ${}^{\mathrm{14}}\mathrm{C}$ values attributed to in situ cosmogenic production (Fang et al., 2021; Hoffmann, 2016). A laboratory irradiation experiment of glacier ice with an artificial neutron flux showed that 11 %–25 % of produced ${}^{\mathrm{14}}\mathrm{C}$ form organic compounds (Hoffmann, 2016). Earlier work involving the irradiation of ice samples to produce ${}^{\mathrm{14}}\mathrm{C}$ (e.g., Roessler et al., 1984) also found that organics accounted for a minor fraction of the total ${}^{\mathrm{14}}\mathrm{C}$. However, we are not aware of any existing studies that specifically investigated production of ${}^{\mathrm{14}}\mathrm{C}$-bearing organic materials in ice from muons. Measuring ${}^{\mathrm{14}}\mathrm{C}$ in organic compounds is unfortunately beyond the scope of this study, as it requires an entirely different analytical setup.

Another possible explanation for the disagreement is that the 2D ice flow model (and thus the time-dependent exposure history of the ice parcels) might be inaccurate. In the following, we conducted sensitivity analyses to combine both uncertainties by using the $+\mathrm{2.4}\mathit{\sigma }$ maximum ablation rate scenario (which corresponds to the deepest physically possible ice trajectory; Fig. S8) and scaling our total ${}^{\mathrm{14}}\mathrm{C}$ upward by 25 % to account for the in situ produced ${}^{\mathrm{14}}\mathrm{C}$ in organics (red dots in Fig. 6d). First, we kept $f_{\mathrm{tot}}$ as a constant ($f_{\mathrm{tot}}=\mathrm{0.021}$, which is the minimum $f_{\mathrm{tot}}$ from the reported uncertainty in Heisinger et al., 2002b) and tuned ${\mathit{\sigma }}_{\mathrm{0}}$ to fit the measurements (dashed red line in Fig. 6d) under the high ablation rate scenario. We find that the best-fit ${\mathit{\sigma }}_{\mathrm{0}}$ (and in situ ${}^{\mathrm{14}}\mathrm{C}$ production from fast muon) is zero. The modeled total ${}^{\mathrm{14}}\mathrm{C}$ under this scenario underestimates the total ${}^{\mathrm{14}}\mathrm{C}$ at lower depths where production from fast muon dominates ($>\mathrm{20}$ m depth) and overestimates the total ${}^{\mathrm{14}}\mathrm{C}$ at depths where production from negative muon capture dominates ($<\mathrm{20}$ m depth). This means that, even under these extreme scenarios, the ${}^{\mathrm{14}}\mathrm{C}$ production rate from negative muon capture has to be lower than the lower-bound estimate of Heisinger et al. (2002b), and some production from fast muons is needed to compensate for the lower production rate from negative muon capture to improve the fit. We then repeated our grid search approach (Sect. 3.2) to find the best-fit $f_{\mathrm{tot}}$ and ${\mathit{\sigma }}_{\mathrm{0}}$ that correspond to the 25 % higher total ${}^{\mathrm{14}}\mathrm{C}$ and maximum ablation rate scenario (solid red line in Fig. 6d). The best-fit $f_{\mathrm{tot}}$ is 0.0055, which is 22 % of the Heisinger et al., 2002b, value) and the best-fit ${\mathit{\sigma }}_{\mathrm{0}}$ is 0.0040 mb (which is 46 % of the Heisinger et al., 2002a, value). We note that these values are within uncertainties in our original best-fit $f_{\mathrm{tot}}$ and ${\mathit{\sigma }}_{\mathrm{0}}$ (Table 2) and still cannot be reconciled with values from Heisinger et al. (2002a, b). We conclude that additional uncertainties from the ice flow history and ${}^{\mathrm{14}}\mathrm{C}$ contribution from organics likely cannot reconcile the difference between the negative muon capture ${}^{\mathrm{14}}\mathrm{C}$ production rate inferred by our data and that of Heisinger et al. (2002b).

In their experimental determination of the ${}^{\mathrm{14}}\mathrm{C}$ production rate by fast muons, Heisinger et al. (2002a) used a single muon energy of 190 GeV ($\mathit{\sigma }\left(E\right)$). The reference nuclear reaction cross section at 1 GeV (${\mathit{\sigma }}_{\mathrm{0}}$) was then scaled using the following equation: 12$$\mathit{\sigma }\left(E\right)={\mathit{\sigma }}_{\mathrm{0}}{\stackrel{\mathrm{‾}}{E}}^{\mathit{\alpha }},$$ where $\mathit{\alpha }$ is a power factor that describes the energy dependence of the cross section (unitless). However, the mean muon energy ($\stackrel{\mathrm{‾}}{E}$) of 190 GeV used by Heisinger et al. (2002a) and the muon flux intensity were much higher than those expected in the first few hundred meters of ice in natural settings (for the top 200 m of Taylor Glacier, $\stackrel{\mathrm{‾}}{E}=\mathrm{32}$ GeV; Fig. S10). It may be possible that the power factor $\mathit{\alpha }$ of 0.75 is incorrect. Balco (2017) has tried fitting the calibration data with $\mathit{\alpha }=\mathrm{1}$ (which simplifies Eq. 12 into a linear relationship between $\mathit{\sigma }$ and $E$). Following Balco (2017), we conducted a sensitivity analysis, fitting our data with similar methods to those described above but with $\mathit{\alpha }=\mathrm{1}$. Increasing $\mathit{\alpha }$ from 0.75 to 1, while keeping ${\mathit{\sigma }}_{\mathrm{0}}$ constant, reduces the overall ${}^{\mathrm{14}}\mathrm{C}$ production rate from fast muons (Fig. S13). To compensate for the lower production rate from fast muons, the new best-fit ${\mathit{\sigma }}_{\mathrm{0}}$ is now 0.0032 mb, which is 29.8 % higher than the best-fit ${\mathit{\sigma }}_{\mathrm{0}}$ when $\mathit{\alpha }=\mathrm{0.75}$ and 36.6 % that of Heisinger et al. (2002b) value. The best-fit $f_{\mathrm{tot}}$ (overall probability of ${}^{\mathrm{14}}\mathrm{C}$ production from negative muon capture) becomes 0.0051, which is 17.2 % higher than the best-fit $f_{\mathrm{tot}}$ when $\mathit{\alpha }=\mathrm{0.75}$ and 20.4 % that of Heisinger et al. (2002a) value (Fig. 7a). These values are within the uncertainties in the $f_{\mathrm{tot}}$ and ${\mathit{\sigma }}_{\mathrm{0}}$ derived with $\mathit{\alpha }=\mathrm{0.75}$ that we presented (Table 2; Fig. 7a) and still cannot be reconciled with original Heisinger et al. (2002a, b) values.

5 Conclusions

This study presents ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements in ablating ice obtained via a new ice sublimation technique, combined with ${}^{\mathrm{14}}\mathrm{CO}$ and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$ measurements obtained from a well-established large-volume melt–extraction method to estimate the species-specific and total in situ muogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates in ice. Under the assumption that the majority of in situ ${}^{\mathrm{14}}\mathrm{C}$ in ice exists as ${}^{\mathrm{14}}\mathrm{CO}$, ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$, and ${}^{\mathrm{14}}{\mathrm{CH}}_{\mathrm{4}}$, we estimated lower muogenic in situ ${}^{\mathrm{14}}\mathrm{C}$ production rates in ice (by a factor of 5.7 (3.6–13.9) and 3.7 (2.0–11.9) with 95 % CI for negative muon capture and fast muon interactions, respectively) compared to values inferred from laboratory irradiation experiments (Heisinger et al., 2002a, b) and measurements in quartz (Lupker et al., 2015; Balco 2017). Prior ice core measurements from Scharffenbergbotnen (Van der Kemp et al., 2002) obtained with an independent technique also appear to be consistent with these lower muogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates in ice. This comparison with the Van der Kemp et al. (2002) results also suggests that a dry extraction technique appears to release essentially all in situ ${}^{\mathrm{14}}\mathrm{C}$ in bubbly (non-clathrated) ice.

At present, there does not appear to be a way to reconcile our Taylor Glacier ice core results and the independent ice core measurements from Scharffenbergbotnen (Van der Kemp et al., 2002) with the muogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates determined in quartz (Heisinger et al., 2002a, b; Lupker et al., 2015). This is a problem that needs further investigation, and we recommend that future studies address this via laboratory muon irradiation experiments involving both ice and quartz targets and studies that include a quantification of the organic fraction of muogenic ${}^{\mathrm{14}}\mathrm{C}$ in ice.

Finally, the constraints on muogenic ${}^{\mathrm{14}}\mathrm{C}$ production rates in ice and the partitioning between the in situ produced ${}^{\mathrm{14}}\mathrm{C}$-bearing gas species provided by this study will allow for future measurements of ${}^{\mathrm{14}}\mathrm{C}$-containing gases in other ice cores to be used for several applications, including using ${}^{\mathrm{14}}{\mathrm{CO}}_{\mathrm{2}}$ measurements for the absolute dating of the bubbles in ice cores (Andree et al., 1984; Van De Wal et al., 1994) and using ${}^{\mathrm{14}}\mathrm{CO}$ measurements to either constrain the past oxidative capacity of the atmosphere (Brenninkmeijer et al., 1992; Petrenko et al., 2021) or reconstruct the past cosmic ray flux (BenZvi et al., 2019).

Data availability

Data from this work will be made available through the U.S. Antarctic Program Data Center (https://www.usap-dc.org/view/project/p0000283, Dyonisius and Petrenko, 2023).

The supplement related to this article is available online at: https://doi.org/10.5194/tc-17-843-2023-supplement.

Author contributions

VVP, EJB, and JPS designed the study. MND and VVP conducted field logistical preparations. MND, VVP, PN, AMS, JAM, SAS, HR, BB, EJB, and JRM conducted the field sampling, on-site sample cutting, and processing. MND, VVP, AMS, and PN extracted the large air samples using the on-site large-volume melter. MND extracted the ${\mathrm{CH}}_{\mathrm{4}}$ and CO from large air samples. BH, VVP, MND, and PFP developed and tested the sublimation system with input from JS. MND extracted the CO${}_{\mathrm{2}}$ using the newly built sublimation system, with assistance from BH and PFP. QH and BY graphitized the ${}^{\mathrm{14}}\mathrm{C}$ samples. AMS conducted the ${}^{\mathrm{14}}\mathrm{C}$ measurements. CA and JAM conducted the ${\mathrm{CH}}_{\mathrm{4}}$ CFA measurements under the supervision of EJB. CA developed the age scale under the supervision of EJB. SEM and IV made the ${\mathit{\delta }}^{\mathrm{13}}$C-${\mathrm{CH}}_{\mathrm{4}}$ stable isotopes measurements. JPS did the $\mathrm{Xe}/\mathrm{Kr}$, $\mathrm{Kr}/{\mathrm{N}}_{\mathrm{2}}$, and $\mathrm{Xe}/{\mathrm{N}}_{\mathrm{2}}$ measurements. RB made the ${\mathit{\delta }}^{\mathrm{15}}$N${}_{\mathrm{2}}$, ${\mathit{\delta }}^{\mathrm{18}}$O${}_{\mathrm{atm}}$, ${}^{\mathrm{40}}\mathrm{Ar}{/}^{\mathrm{36}}\mathrm{Ar}$, ${\mathrm{O}}_{\mathrm{2}}/{\mathrm{N}}_{\mathrm{2}}$, and $\mathrm{Ar}/{\mathrm{N}}_{\mathrm{2}}$ measurements. CH did the [${\mathrm{CH}}_{\mathrm{4}}$] and halogenated trace gas measurements under the supervision of RFW. MK made the discrete [${\mathrm{CH}}_{\mathrm{4}}$] mole fraction and total air content measurements. IV did the ${\mathit{\delta }}^{\mathrm{13}}$CO measurement for the CO dilution gas. CB developed the ice flow model. MND developed the ${}^{\mathrm{14}}\mathrm{C}$ production model, with input from CB and VVP. MND, BH, and VVP analyzed the results and wrote the paper, with input from all co-authors.

Competing interests

The contact author has declared that none of the authors has any competing interests.

Disclaimer

Publisher’s note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Acknowledgements

We thank Mike Jayred of the U.S. Ice Drilling Program (IDP), for ice drilling, camp manager Kathy Schroeder, for assistance in the field, and the United States Antarctic Program (USAP), for logistical support. We thank Emily Mesiti, for her assistance in the Rochester ice core lab during the sublimation extraction campaign. This work has been supported by the U.S. NSF awards (grant nos. PLR-1245659 for VVP, PLR-1245821 for EJB, and PLR-1246148 for JPS) and the Packard Fellowship for Science and Engineering (VVP).

Financial support

This research has been supported by the Office of Polar Programs (grant nos. PLR-1245659, PLR-1245821, and PLR-1246148).

Review statement

This paper was edited by Benjamin Smith and reviewed by Greg Balco, Nathaniel A. Lifton, and one anonymous referee.