Last ice sheet recession and landscape emergence above sea level in east-central Sweden, evaluated using in situ cosmogenic ¹⁴C from quartz

In situ cosmogenic ¹⁴C (in situ ¹⁴C) in quartz provides a recently developed tool to date exposure of bedrock surfaces of up to ∼ 25 000 years. From outcrops located in east-central Sweden, we tested the accuracy of in situ ¹⁴C dating against (i) a relative sea level (RSL) curve constructed from radiocarbon dating of organic material in isolation basins and (ii) the timing of local deglaciation constructed from a clay varve chronology complemented with traditional radiocarbon dating. Five samples of granitoid bedrock were taken along an elevation transect extending southwestwards from the coast of the Baltic Sea near Forsmark. Because these samples derive from bedrock outcrops positioned below the highest postglacial shoreline, they target the timing of progressive landscape emergence above sea level. In contrast, in situ ¹⁴C concentrations in an additional five samples taken from granitoid outcrops above the highest postglacial shoreline, located 100 km west of Forsmark, should reflect local deglaciation ages. The 10 in situ ¹⁴C measurements provide robust age constraints that, within uncertainties, compare favourably with the RSL curve and the local deglaciation chronology. These data demonstrate the utility of in situ ¹⁴C to accurately date ice sheet deglaciation, and durations of postglacial exposure, in regions where cosmogenic ¹⁰Be and ²⁶Al routinely return complex exposure results.

1 Introduction

The pacing of the retreat of ice sheets in North America and Eurasia, since their maximum expansion during the last glaciation, remains an active research field (e.g. Hughes et al., 2016; Stroeven et al., 2016; Patton et al., 2017; Dalton et al., 2020, 2023). Understanding the triggers and processes causing the demise of these ephemeral ice sheets yields the best blueprint for understanding the future behaviour of the Greenland and Antarctic ice sheets in a warming climate. Coupling the behaviour of deglaciating ice sheets over the course of the Late Glacial and Early Holocene to increasingly precise climate reconstructions, including climatic events, requires increased precision in ice sheet reconstructions (e.g. Bradwell et al., 2021). Precision can be enhanced through coupling geomorphological mappings of ice sheet margins (such as moraines, grounding zone wedges, lateral meltwater channels, and ice-dammed lake shorelines and spillways) with numerical field constraints from a diverse array of dating techniques (e.g. Stroeven et al., 2016; Bradwell et al., 2021; Regnéll et al., 2023).

Ice sheet reconstructions, especially in North America, have become highly detailed through radiocarbon dating (Dyke et al., 2002; Dalton et al., 2020). With the advance of offshore imaging of glacial geomorphology (Greenwood et al., 2017, 2021; Bradwell et al., 2021), radiocarbon dating has received a renewed upswing in recent years (e.g. Dalton et al., 2020; Bradwell et al., 2021). However, large landscape areas lack radiocarbon age constraints on ice sheet retreat because of an absence of datable organic material. Fortunately, optically stimulated luminescence ages on buried sand layers (e.g. Alexanderson et al., 2022) and cosmogenic nuclide-apparent exposure ages on exposed bedrock and erratics have narrowed some of the gaps (e.g. Hughes et al., 2016; Stroeven et al., 2016; Dalton et al., 2023). In studies using cosmogenic nuclides, an “apparent” exposure age is derived from a simple calculation from the nuclide concentration under consideration (Lal, 1991; Gosse and Phillips, 2001). Correctly interpreting the exposure age relies on modelling that considers geological factors that can reduce the nuclide concentration relative to the time since initial subaerial exposure (such as erosion and burial by glacial ice, water, snow, and/or soil; Gosse and Phillips, 2001; Schildgen et al., 2005; Ivy-Ochs and Kober, 2008). Exposure dating is the only technique available in regions where ice sheet erosion has left the surface bare or covered by a thin drape of till. Kleman et al. (2008) show that for Fennoscandia, these conditions are widespread in coastal regions where ice accelerated towards its streaming sectors and where wave wash during glacial rebound further thinned or removed pre-existing sediment covers.

Coastal sectors in formerly glaciated regions provide sites important to the study of paleoglaciology. They offer an abundance of bedrock exposures from which patterns and processes of subglacial erosion can be studied through cosmogenic nuclide exposure dating (e.g. Hall et al., 2020). Also, because of the interplay with postglacial sea levels, coastal areas yield data on glacio-isostatic rebound that are critical to the geodynamic modelling of Earth rheology and thicknesses of former ice sheets (e.g. Lambeck et al. (1998, 2010) and Patton et al. (2017) for Fennoscandian examples). Geodynamic models require validation against measurements of vertical crustal motion (Steffen and Wu, 2011), such as those provided by recent global positioning system (GPS) measurements (e.g. Lidberg et al., 2010) and postglacial records of crustal rebound afforded by relative sea level (RSL) curves (e.g. Påsse and Andersson, 2005). The construction of RSL curves, detailing the history of land surface emergence from sea level, is traditionally done either using sediments accumulated in isolation basins at different elevations above sea level or by dating uplifted gravel beach ridges. Typically, isolation basins, and their sediments, show a progression from marine to brackish and finally to freshwater environments as they are uplifted through tidal levels (Long et al., 2011). Histories of land uplift above sea level are documented using microfossil and macrofossil analyses of isolation basin sediments and radiocarbon dating on macrofossils (Romundset et al., 2011). Uplifted beach ridges can be radiocarbon-dated from a variety of materials (Blake, 1993) but most confidently from driftwood, whalebone, and shells (e.g. Dyke et al., 1992). Gravel beach ridges have also been investigated using OSL (optically stimulated luminescence) and ¹⁰Be exposure dating even though, other than the highest beach ridge, they may be prone to clast reworking (Briner et al., 2006; Simkins et al., 2013; Bierman et al., 2018). A distinct advantage of constructing RSL curves using cosmogenic nuclides is that land surface emergence above sea level may be additionally dated from boulders (Briner et al., 2006) or bedrock (Bierman et al., 2018).

The potential for cosmogenic surface exposure dating of the last ice sheet retreat in recently glaciated low-relief cratonic landscapes would seemingly be high because of the frequent outcropping of glacially sculpted quartz-bearing crystalline bedrock. However, either the ice sheet may have been non-erosive, or the erosion was insufficiently deep to remove all the cosmogenic nuclide inventory from previous exposure periods. Apparent ages are therefore often older than indicated by radiocarbon dating (Heyman et al., 2011; Stroeven et al., 2016) because they include a component of nuclide inheritance. Apparent ages younger than indicated by radiocarbon dating can also occur if sampled rock surfaces have been shielded, for example by sediments, following deglaciation. Concentrations of ¹⁰Be and ²⁶Al, in either bedrock or erratic boulders, often reflect complex exposure histories rather than simple deglacial exposure durations (Heyman et al., 2011; Stroeven et al., 2016).

In this study we use ¹⁴C produced in situ in quartz-bearing bedrock (in situ ¹⁴C) because it potentially circumvents an overt reliance on the need for deep erosion (> 3 m) to remove the inherited signal from previous exposure periods (Gosse and Phillips, 2001). Because of its short half-life of 5700 ± 30 years, inherited in situ ¹⁴C will decay if ice sheet burial at investigated sites during the last glacial phase (marine isotope stage 2; MIS2) exceeded 25–30 kyr, which is ca. five half-lives (Briner et al., 2014).

Some studies assessing changes in glacier and ice sheet extents over Late Glacial–Holocene timescales have used in situ ¹⁴C (Miller et al., 2006; Fogwill et al., 2014; Hippe et al., 2014; Schweinsberg et al., 2018; Pendleton et al., 2019; Young et al., 2021; Schimmelpfennig et al., 2022). In these studies, in situ ¹⁴C has been applied with other nuclides with longer half-lives, in particular ¹⁰Be, to unravel complex histories of glacier advance and retreat (e.g. Goehring et al., 2011) and spatial patterns in glacial erosion in mountainous terrain (e.g. Steinemann et al., 2021). Extensive regions formerly covered by ice sheets are characterized by low-relief and low-elevation terrain. The effectiveness of in situ ¹⁴C in dating ice sheet retreat in these non-alpine settings and in quantifying shoreline displacement from bedrock samples has not been previously assessed. The aim of this study is therefore to validate the use of ¹⁴C formed in situ in bedrock as a reliable chronometer by evaluating its performance in duplicating (i) a previously established Holocene RSL curve based on radiocarbon dating (Hedenström and Risberg, 2003; SKB, 2020) and (ii) the timing of deglaciation above the highest (postglacial) shoreline in nearby east-central Sweden according to reconstructions of deglaciation of the last ice sheet (Hughes et al., 2016; Stroeven et al., 2016).

2 Study area

Our study is focused on east-central Sweden. This region includes low elevation and low relief terrain in Forsmark–Uppland and higher elevation and higher relief terrain in adjacent Dalarna–Gävleborg (Fig. 1). This region was selected because Forsmark is the location of a planned geological repository for spent nuclear fuel (e.g. SKB, 2022). As such, this region has been intensively studied and has a wealth of geologic data relevant to our study. This wealth includes in-depth analyses of bedrock and environmental properties, including influences of glacial and postglacial processes (e.g. Lönnqvist and Hökmark, 2013; Hall et al., 2019; Moon et al., 2020; SKB, 2020).

Figure 1 Sample locations for in situ ¹⁴C dating in (a) Dalarna–Gävleborg and (b) Forsmark–Uppland. The five Dalarna–Gävleborg sample sites are located on what were islands above the highest postglacial shoreline (shown), whereas the five sample sites from Forsmark–Uppland are located below the highest shoreline (not shown because the entire area was submerged). See inset maps for locations of panels (a) and (b) and for the 10.7 and 10.8 ka retreat isochrones (blue) from Stroeven et al. (2016) and 11 ka (most-credible, minimum, and maximum) retreat isochrones (red) from Hughes et al. (2016). The rectangle in panel (b) approximately indicates the site selected for the planned geological repository for spent nuclear fuel at Forsmark. DEM with 2 m resolution, from lidar data (Lantmäteriet, 2020).

From spatio-temporal ice sheet reconstructions by Kleman et al. (2008), the study area was glaciated 16–20 times for a total duration of ca. 330 kyr over the past 1 Myr. The last deglaciation of the study area is well-constrained by two recent reconstructions that differ in their approach (Hughes et al., 2016; Stroeven et al., 2016). The Hughes et al. (2016) reconstruction relies primarily upon chronological constraints supplied from radiocarbon, thermal luminescence, optically stimulated luminescence (OSL), infrared stimulated luminescence, electron spin resonance, terrestrial cosmogenic nuclide (TCN), and U-series dating. Published landform data, mostly with respect to end moraines and generally accepted correlations of ice margin positions between individual moraines, provide complementary evidence. In contrast, the Stroeven et al. (2016) reconstruction combines geomorphological constraints for ice sheet margin outlines, including ice-marginal depositional landforms and meltwater channels, ice-dammed lakes, eskers, lineations, and striae, with chronological constraints supplied by radiocarbon, varve, OSL, and TCN dating. Whereas Hughes et al. (2016) reconstruct ice sheet retreat every 1 kyr and for every ice margin plot its position as “most credible”, “min”, and “max”, Stroeven et al. (2016) present ice margin positions for every 100 years inside the Younger Dryas standstill position (Stroeven et al., 2015). These marginal positions are temporally and spatially defined by the “Swedish timescale” clay varve record along the Swedish east coast (De Geer, 1935, 1940; Strömberg, 1989, 1994; Brunnberg, 1995; Wohlfarth et al., 1995). From Stroeven et al. (2016), the last deglaciation of the study area occurred at 10.8 ± 0.3 ka, which overlaps the timing of deglaciation of the study area from Hughes et al. (2016), within uncertainty (Fig. 1). The highest postglacial shoreline in east-central Sweden is located at a present elevation of ∼ 200 m a.s.l. in Dalarna–Gävleborg, which is ∼ 100 km west of Forsmark (SGU, 2015). The exposure duration of bedrock above the highest postglacial shoreline represents the time since local deglaciation. Hence, in situ ¹⁴C ages from bedrock above the highest postglacial shoreline should conform to the reconstructed deglaciation age of 10.8 ± 0.3 ka from Stroeven et al. (2016).

Below the highest postglacial shoreline, in the Forsmark–Uppland region, the last deglaciation occurred in a marine environment, and the landscape has progressively emerged above sea level through postglacial isostatic uplift. A RSL curve constructed from radiocarbon dating of basal organic sediments trapped in isolation basins along elevation transects describes the progressive emergence of the Forsmark–Uppland landscape above sea level (Robertsson and Persson, 1989; Risberg, 1999; Bergström, 2001; Hedenström and Risberg, 2003; Berglund, 2005; SKB, 2020). Ages calculated from in situ ¹⁴C from bedrock outcrops along an elevation transect would then mirror the Forsmark RSL curve for their corresponding elevations (but they would be slightly older because of nuclide production through shallow water before emergence).

A potential complication to the accurate exposure age dating of bedrock surfaces using in situ ¹⁴C in east-central Sweden is that the most recent period of ice sheet burial may have been insufficiently long to decay in situ ¹⁴C inventory inherited from prior exposure. Here, the extent of the Fennoscandian Ice Sheet during interstadial MIS3 and the timing of ice advance across the Forsmark region during late MIS3 are crucially important. Kleman et al. (2020) have identified ice-free conditions around Idre (330 km NW, up-ice, of our study area; Fig. 1) between 55 and 35 ka, which implies inundation of our study area by ice after 35 ka. Combined with a well-constrained final deglaciation age of 10.8 ± 0.3 ka (Stroeven et al., 2016), it appears that our study area has most recently (during MIS2) been inundated by glacial ice for, at most, 24 kyr. This inference is in line with results from ice sheet modelling indicating a 22 kyr duration of ice cover at Forsmark during MIS2 (SKB, 2020). Consequently, it is possible that in situ ¹⁴C concentrations may reflect subaerial exposure of bedrock in our study area during MIS3, in addition to Holocene exposure, resulting in an offset towards older ages relative to the RSL curve for Forsmark (Hedenström and Risberg, 2003; SKB, 2020) and the deglaciation chronologies of Hughes et al. (2016) and Stroeven et al. (2016).

3 Methods

3.1 Sampling of bedrock outcrops for in situ ¹⁴C measurements

We used the following sampling strategy to evaluate the accuracy of bedrock exposure ages derived from in situ ¹⁴C against the Forsmark RSL curve and the deglaciation of the last ice sheet in east-central Sweden. A rigorous scheme was applied to ensure that we avoided sampling quartz altered through hydrothermal processes that is likely to occur in major pegmatite intrusions; outcrops located in major deformation zones; and outcrop-scale veins, fractures, and adjacent rock volumes. Consequently, sampling was done on outcrops of metagranitoid from the early-Svecokarelian GDG-GSDG suite that dominates the Bergslagen lithotectonic unit (Stephens and Jansson, 2020). A petrological examination using transmitted light polarization microscopy was applied to thin sections to ascertain that the quartz was unlikely to contain multi-fluid-phase, vapour-phase, or solid-phase inclusions. All samples were collected using an angle grinder, which permits the sampling of hard crystalline bedrock isolated from outcrop edges, fractures, and quartz veins and consistently limits sample thicknesses to 3 cm.

We collected a total of 10 samples for in situ ¹⁴C analyses. Five of these were collected along a SW–NE transect near Forsmark (Fig. 1b). These outcrops were chosen because they span an elevation gradient of 9.4–56.0 m a.s.l., and the exposure ages derived from in situ ¹⁴C can therefore be evaluated against the Forsmark RSL curve. We collected a further five samples from locations above the highest shoreline (Fig. 1a) to determine the age of local deglaciation for comparison with published deglaciation chronologies (Hughes et al., 2016; Stroeven et al., 2016). Sample locations were logged on a 2 m resolution lidar digital elevation model (DEM) displayed in ArcGIS 10 on a tablet computer. A GPS add-in tool in ArcGIS 10 was used to record positional data within a horizontal precision of 2 m. The elevation of each sample location was extracted from the DEM and has a precision of tens of centimetres. The influence of these minor positional uncertainties on our ¹⁴C calculations is trivial, and none of the sample sites are influenced by topographic shielding that could reduce the accumulation of ¹⁴C in bedrock.

Each sampled bedrock outcrop formed a local topographic high, which minimizes the risk of burial by soil and snow (Supplement file S1). Moss mats were present on all sampled outcrops. Although we avoided sampling bedrock that was moss-covered, we cannot be certain that moss mats did not formerly cover the sample sites. Given a compressed thickness of 0.5 cm and an estimated density of 0.7 g cm⁻³, this may have contributed to a shielding of the sampled rock surfaces of 0.35 g cm⁻², which is negligible and is therefore excluded from our age inferences.

Figure 2 Production rate calibration of ¹⁴C in quartz. (a) Reference spallation ¹⁴C production rate calibration based on data from Schimmelpfennig et al. (2012), Young et al. (2014), Lifton et al. (2015), Borchers et al. (2016), and Phillips et al. (2016), corrected per Hippe and Lifton (2014) and compiled in Koester and Lifton (2023). An uncertainty-weighted production rate is calculated for each of the eight sites. Outliers, which are not included in the uncertainty-weighted production rates, are determined based on the requirement that there should be at least three samples yielding a reduced chi-square statistic $\left(X_{\mathrm{R}}^{\mathrm{2}}\right)$ with a p value of at least 0.05 for the assumption that the individual production rates from a site are derived from one normal distribution. For $X_{\mathrm{R}}^{\mathrm{2}}$, but not the uncertainty weighting, we use the largest of the sample-specific production rate uncertainty based on the ¹⁴C concentration uncertainties and 5 % of the sample production rate. This procedure does not punish samples with low measurement uncertainties, which otherwise risk exclusion as outliers. We adopt a global reference spallation ¹⁴C production rate of 12.81 ± 1.25 atoms g⁻¹ yr⁻¹, calculated as the arithmetic mean of the eight site production rates with the uncertainty being based on an uncertainty-weighted deviation of all included single-sample production rates, excluding outliers. (b) Calibration of ¹⁴C production rate from muons based on the data of Lupker et al. (2015). The calibration is based on the method used in the CRONUScalc calculator (Marrero et al., 2016; Phillips et al., 2016). The figure shows the best-fit ¹⁴C concentration profiles produced from spallation, slow muons, and full production. The best fit yields a near-zero production from fast muons (cf. Lupker et al., 2015). The production rate calibration has been carried out using the expage calculator version 202403 (Heyman, 2024) in an iterative way to make the global reference spallation ¹⁴C production rate converge with the production rate from muons.

3.2 Laboratory preparation for accelerator mass spectrometry (AMS)

Samples were physically and chemically processed at the Purdue Rare Isotope Measurement Laboratory (PRIME Lab) at Purdue University, USA. Concentrations of in situ ¹⁴C were determined from purified quartz separates through automated procedures (Lifton et al., 2023). Approximately 5 g of quartz from each sample was added to a degassed LiBO₂ flux in a reusable 90 % Pt$/$10 % Rh sample boat and heated to 500 °C for 1 h in ca. 6.7 kPa of research purity O₂ to remove atmospheric contaminants, which were discarded. The sample was then heated to 1100 °C for 3 h to dissolve the quartz and release the in situ ¹⁴C, again in an atmosphere of ca. 6.7 kPa of research purity O₂ to oxidize any evolved carbon species to CO₂. The CO₂ from the 1100 °C step was then purified, measured quantitatively, and converted to graphite for ¹⁴C AMS measurement at PRIME Lab (Lifton et al., 2023). To test for data reproducibility, sample BG21-002 was randomly selected to undergo laboratory preparation and AMS a second time. Measured concentrations of in situ ¹⁴C are calculated from the measured isotope ratios via AMS following Hippe and Lifton (2014; Table 1).

Figure 3 Apparent ¹⁴C exposure ages for five Forsmark samples from below the highest shoreline (Fig. 1b; Table 2) with 1σ external uncertainties. The expected exposure ages are calculated assuming the RSL curve is correct, the ¹⁴C spallation production rate is correct, partial exposure as the sample approaches the water surface, and full postglacial exposure for the duration above sea level. Hence, the expected exposure age curve is a few hundred years older than the RSL curve. The RSL curve is from SKB (2020), and uncertainties for the 1–6 ka interval are calculated from the original radiocarbon data in Hedenström and Risberg (2003). The RSL uncertainty envelope is also transposed onto the expected exposure age curve.

Figure 4 Exposure ages from samples above the highest shoreline (Fig. 1a; Table 2). The individual samples (filled black circles) display a 1σ internal uncertainty (measurement uncertainty; black lines). For the repeat sample BG21-002, the exposure age is calculated with a weighted mean ¹⁴C concentration using a 2 % uncertainty. The cosmogenic nuclide ages yield a reduced chi-square $\left(X_{\mathrm{R}}^{\mathrm{2}}\right)$ statistic of 1.78 and a p value of 0.13 based on internal uncertainties, which indicates that they are from the same population. The colour gradient for the Stroeven et al. (2016) deglaciation chronology indicates the 0.1–0.5 ka uncertainty range, whereas the uncertainty for the Hughes et al. (2016) chronology reflects the maximum and minimum estimates for deglaciation of the study area, which are unequally distributed around the most credible estimate (orange line).

3.3 Exposure age calculations

The expage calculator version 202403 (Heyman, 2024) is used to calculate apparent exposure ages. It is based on the original version 2 CRONUS calculator (Balco et al., 2008), the LSDn production rate scaling (Lifton et al., 2014), and the CRONUScalc calculator (Marrero et al., 2016), using the geomagnetic framework of Lifton (2016) with the SHA.DIF.14k model for the last 14 kyr. Exposure ages are calculated using resulting time-varying ¹⁴C production rates accounting for decay and interpolated to match the measured ¹⁴C concentration. The production rate from muons is calibrated against the Leymon High core ¹⁴C data of Lupker et al. (2015), and the production rate from spallation is calibrated against updated global ¹⁴C production rate calibration data (Schimmelpfennig et al., 2012; Young et al., 2014; Lifton et al., 2015; Borchers et al., 2016; Phillips et al., 2016; Koester and Lifton, 2023, 2024). This calibration is done iteratively for spallation and muons to reach convergence using the expage production rate calibration methods (Fig. 2; Heyman, 2024).

Exposure age calculations along the Forsmark–Uppland transect account for ¹⁴C production during emergence through shallow water. Burial of sampled surfaces by snow is excluded from the age calculations for all sample sites because we know neither how snow burial depths and durations vary between sites nor how they vary through time. The effect of snow burial would be to slightly decrease cosmogenic nuclide production in the underlying rock surface (Schildgen et al., 2005), and we have minimized this effect through our sampling strategy.

Table 1 In situ ¹⁴C sample measurement details.

¹ Purdue Carbon Extraction and Graphitization System. ² PRIME Lab ID. ³ Mass graphitized for AMS analysis after a small aliquot (ca. 9 µg C) was taken for an offline stable C isotopic analysis. ⁴ VPDB stands for Vienna Peedee Belemnite. ⁵ Measured relative to OX-2 standard. ⁶ Corrected for mass-dependent graphitization blank (based on AMS split mass C) and stable C composition. ⁷ Sample values calculated using diluted mass C and corrected for mean procedural blank (all blanks).

Figure 5 Modelled in situ ¹⁴C concentration evolution over the last 80 kyr in the five samples (BG21-001–BG21-005) from above the highest shoreline. The ¹⁴C development is modelled assuming no glacial or interglacial erosion, continuous exposure to cosmic rays during ice-free periods, and full shielding from cosmic rays (no ¹⁴C production) during periods with ice cover. The points just left of the plots display the measured ¹⁴C concentrations for the six sample measurements (Table 1). (a) Scenario with short periods of MIS4 and MIS2 ice cover from 66 to 60 ka and from 28 ka to deglaciation around 10.7 ka. (b) Scenario with longer periods of MIS4 and MIS2 ice cover from 70 to 57 ka and from 35 ka to the deglaciation around 10.7 ka. Due to the rapid decay of ¹⁴C (half-life of 5700 ± 30 years), both scenarios yield similar end-point concentrations of ¹⁴C that overlap, within uncertainties, the measured sample concentrations.

4 Results

Analytical results for in situ ¹⁴C samples and procedural blanks are presented in Table 1. The mean and standard deviation are used to correct measured ¹⁴C sample inventories (Table 1) because procedural blanks are well-constrained during the analytical time frame. Inferred ages for the five in situ ¹⁴C samples from the Forsmark–Uppland transect (i.e. below the highest postglacial shoreline) are shown relative to the Holocene RSL curve for Forsmark and the expected in situ ¹⁴C exposure age curve considering subaqueous cosmogenic nuclide production (Fig. 3; Table 2). Exposure age uncertainties are large with internal uncertainties (measurement uncertainties; Balco et al., 2008) of 5 %–9 % and external uncertainties of 13 %–25 % (also including production rate uncertainties, which are high relative to ¹⁰Be; Borchers et al., 2016; Phillips et al., 2016). Apparent exposure ages increase consistently with elevation and match expected ages within uncertainty. The two highest samples have near-identical apparent exposure ages and elevations. However, these samples provide independent ages because they are horizontally separated by 624 m (Fig. 1b). There is good agreement between ages inferred from these in situ ¹⁴C data and the RSL curve constructed from organic radiocarbon dating of isolation events (Hedenström and Risberg, 2003; SKB, 2020).

Apparent exposure ages for the five in situ ¹⁴C samples located above the highest shoreline in Dalarna and Gävleborg (Fig. 1a) are shown in Fig. 4 and Table 2. The weighted mean age from all five samples is 11.9 ± 1.5 ka. These data display a $X_{\mathrm{R}}^{\mathrm{2}}$ of 1.78 and a p value of 0.13 based on 1σ internal uncertainties, which does not support a rejection of the hypothesis that the apparent exposure ages represent the same population. In addition to the samples being from the same population, the exposure ages are consistent, within uncertainty, with the expected deglaciation age of 10.8 ± 0.3 ka (Stroeven et al., 2016). Replicate measurements of sample BG21-002 closely agree, and an age based on a weighted mean ¹⁴C concentration is shown in Fig. 4.

Table 2 Apparent in situ ¹⁴C ages from bedrock quartz located in Dalarna–Gävleborg and Forsmark–Uppland.

¹ All samples have a thickness of 3 cm, a density of 2.7 g cm⁻³, and a shielding factor of 1. Zero erosion is assumed. ² The ¹⁴C age and 1σ external uncertainty (1σ internal uncertainty).

5 Discussion

The in situ ¹⁴C bedrock exposure ages from the Forsmark–Uppland transect (i.e. below the highest postglacial shoreline) consistently increase with elevation and overlap the expected exposure age curve, within uncertainty (Fig. 3). This study adds to the few previous applications of cosmogenic nuclides to defining postglacial landscape emergence above sea level (Briner et al., 2006; Bierman et al., 2018). Briner et al. (2006) present good (visual) congruence with a record of shoreline emergence built from radiocarbon-dated driftwood and fauna by Dyke et al. (1992) using ¹⁰Be measurements in boulders on beaches derived from wave-washed till. Their study also mentions that building a relative sea level curve from pebbles, cobbles, and plucked bedrock suffered from inheritance problems, an experience shared by Matmon et al. (2003) while attempting the dating of chert on beach ridges in southern Israel and heeded by Bierman et al. (2018). Bierman et al. (2018) successfully dated landscape emergence on Greenland using ¹⁰Be across a range of settings, including bedrock below the highest shoreline, cobbles from beach ridges at the highest shoreline, and boulders and bedrock above the highest shoreline. They note that success hinges on the requirement of warm-based ice and deep glacial erosion in exposing bedrock devoid of an inherited cosmogenic nuclide inventory. In many regions, however, including east-central Sweden and more widely in Fennoscandia, these conditions are not met either because of cold-based conditions (Patton et al., 2016; Stroeven et al., 2016) or because of weakly erosive warm-based ice, such as at Forsmark (Hall et al., 2019; SKB, 2020), during all or much of glacial time. Cosmogenic nuclide inheritance is therefore a part of the landscape fabric. Bierman et al. (2018) advocate the use of in situ ¹⁴C as a methodology to circumvent inheritance problems. Our study is the first to follow up on that suggestion and shows, convincingly, that using in situ ¹⁴C can extend the study of landscape rebound to regions where ice sheet erosion was insufficiently deep to allow for the application of long-lived nuclides.

Five bedrock samples from above the highest postglacial shoreline are well-clustered, and the weighted mean age (and full uncertainty) of 11.9 ± 1.5 ka overlaps with the predicted deglaciation age of 10.8 ± 0.3 ka (Fig. 4; Hughes et al., 2016; Stroeven et al., 2016). Because derived exposure ages overlap with the predicted deglaciation age, we further infer that the in situ ¹⁴C samples, including those located below the highest postglacial shoreline, within uncertainty, lack significant inheritance from previous exposure. Model scenarios of in situ ¹⁴C concentration evolution over varying durations of MIS2 and MIS4 ice cover are consistent with minor inheritance, even with short periods of ice coverage and no glacial or interglacial erosion (Fig. 5). Even if the last ice sheet had advanced over the region as late as 28 ka, there would only be a negligible inventory of inherited ¹⁴C atoms produced prior to the MIS2 ice advance.

Our in situ ¹⁴C data from above the highest (postglacial) shoreline demonstrate their potential for constraining the deglaciation chronology of former ice sheets. This is especially true for regions with thin till drapes, abundant bedrock exposures, and sparse moraines outlining successive retreat stages. In Fennoscandia, thin tills commonly occur (cf. Kleman et al., 2008), and ice sheet retreat appears to have proceeded uninterrupted inside the Younger Dryas moraine belt (apart from the Central Finland Ice Marginal Formation; e.g. Rainio et al., 1986; Stroeven et al., 2016). Even though the post-Younger Dryas deglaciation of east-central Sweden is well-constrained by a clay varve chronology below the highest postglacial shoreline (Strömberg, 1989), there are vast areas above the highest shoreline that remain poorly constrained by data (Stroeven et al., 2016). In addition to a lack of datable deglacial landforms, this is attributable to glacial erosion of bedrock having frequently been insufficient to remove inventories of long half-life ¹⁰Be and ²⁶Al (Patton et al., 2022), thereby leaving nuclides inherited from exposure prior to the last glaciation (Heyman et al., 2011; Stroeven et al., 2016). Because of the short ¹⁴C half-life and an improved sampling methodology, in situ ¹⁴C may now be a prime candidate nuclide to be included in last deglaciation studies on glaciated cratons, such as the dating of boulders deposited along glacial flowlines. This is a technique practised successfully using ¹⁰Be (Margold et al., 2019; Norris et al., 2022).

6 Conclusions

A total of 10 in situ ¹⁴C measurements of bedrock are consistent with a RSL curve for Forsmark derived from organic radiocarbon dating of basal sediments in isolation basins and the Fennoscandian Ice Sheet deglaciation chronologies from Stroeven et al. (2016) and Hughes et al. (2016). This study introduces the use of in situ ¹⁴C in Fennoscandian Ice Sheet paleoglaciology and outlines a promise of its use as a basis for supporting future shoreline displacement studies and for tracking the deglaciation in areas that lack datable organic material and where ¹⁰Be and ²⁶Al routinely return complex exposure results.