Progression of the surge in the Negribreen Glacier System from two years of ICESat-2 measurements

Abstract

The state-of-the-art measurement capabilities of ICESat-2 allow high spatiotemporal resolution of complex ice-dynamic processes that occur during a surge. Detailed and precise mapping of height changes on surging glaciers has previously escaped observations from space due to the limited resolution of space-borne altimeter data and the surface characteristics of glaciers during surge, such as heavy crevassing. This makes geophysical interpretation of deformation and assessment of mass transfer difficult. In this paper, we present an approach that facilitates analysis of the evolution of geophysical processes during a surge, including height changes, crevassing, mass transfer and roughness. We utilize all data from two years of ICESat-2 observations collected during the mature phase of the Negribreen Glacier System surge in 2019 and 2020. The progression of Negribreen's surge has resulted in large-scale elevation changes and wide-spread crevassing, making it an ideal case study to demonstrate ICESat-2 measurement capabilities, which are maximized when coupled with the Density Dimension Algorithm for ice surfaces (DDA-ice). Results show expansion of the surge in upper Negribreen which demonstrates the ability of ICESat-2/DDA-ice to measure a rapidly changing surging glacier and provide the best estimates for cryospheric changes and their contributions to sea-level rise.

1. Introduction

1.1 Glacier surging: relevance, complexity and observation from space

Processes of ice-dynamic instability, such as glacial acceleration, have been identified by the Intergovernmental Panel of Climate Change (IPCC, AR6, Chen and others, 2021) as a source of ‘deep uncertainty’ in our understanding of ice sheet and glacier evolution and their contribution to future sea-level rise. Glacier surging is a type of glacial acceleration (Clarke, 1987; Truffer and Echelmeyer, 2003), and is incompletely understood because surges are relatively rare events with dramatic ice deformation that occurs on short timescales, which results in a paucity of comprehensive observations relative to other acceleration types. While more extensive surge observations have been conducted with the advance of satellite technology (e.g. Dunse and others, 2015; Trantow and Herzfeld, 2016; Goerlich and others, 2020; Paul and others, 2021; Kääb and others, 2023; Trantow and Herzfeld, 2024), few of these are of a resolution to detect surface-structure signatures of dynamic change, especially with respect to crevassing, and none so far have used ICESat-2. Here we employ the unique capabilities of ICESat-2 to retrieve surface heights of crevassed ice with the detail needed to infer dynamic and kinematic information of a surging glacier from morphological cross sections.

A surge-type glacier will cycle quasi-periodically between a long quiescent phase of normal flow and a short surge phase when flow speeds rapidly accelerate by a factor of 10–200. Negribreen, a large glacier in eastern Spitsbergen, Svalbard, began to surge in 2016 for the first time in over 80 years (Lefauconnier and Hagen, 1991; Strozzi and others, 2017; Herzfeld and others, 2021). Negribreen continues to surge through 2024, though it is gradually decelerating from its peak velocities in 2017 (Haga and others, 2020; Herzfeld and others, 2022). A feed-back system of acceleration, crevassing and ice-ocean interaction has already led to accelerated mass loss of the Negribreen Glacier System. There is a possibility that the Negribreen surge may trigger pervasive mass loss and a potential disintegration of the entire glacier system, similar to the collapse of the Vavilov Ice Cap in the Russian Arctic (Willis and others, 2018).

Surges are limited geographically to particular regions of the cryosphere, clustered mostly in Alaska/Yukon, Svalbard, high mountain Asia and around the periphery of the major ice sheets (Dolgushin and Osipova, 1975; Kamb and others, 1985; Herzfeld and Mayer, 1997; Björnsson and others, 2003; Jiskoot and others, 2003; Jiskoot, 2011; Flowers and others, 2011; Sevestre and others, 2015; Bhambri and others, 2017; Trantow and Herzfeld, 2018; Kochtitzky and others, 2020; Vale and others, 2021; Banerjee and others, 2022; Guillet and others, 2022; Yao and others, 2023). While new satellite technologies have allowed detection of over 100 surge-type events from 2017 to 2022 (Kääb and others, 2023), very few of these events have received dedicated studies that map and document individual surges, making investigation of surge mechanisms and processes a data-starved problem.

Geophysical processes that are characteristic of a surge include rapid acceleration, crevassing, mass transfer within the glacier, advance of the ice front and changes in the internal hydrological system of the glacier (Meier and Post, 1969; Kamb, 1987; Harrison and Post, 2003; Truffer and others, 2021). Altimetry is capable of measuring these first four surge characteristics in particular. However, observations of surging glaciers from space typically use radar or optical imagery (e.g. Guan and others, 2022; Wuite and others, 2022; Kääb and others, 2023; Liu and others, 2024; Main and others, 2024), which usually focus on velocity or elevation-change analysis yet lack the spatial resolution to resolve crevasses, and field measurements (e.g. Kamb and others, 1985; Björnsson and others, 2003; Herzfeld and others, 2013b, 2022; Lovell and Fleming, 2023), which are relatively rare and often lack spatial and/or temporal coverage.

High-resolution optical imagery from Maxar's Worldview satellites has allowed surge analysis at the meter to submeter resolution (Herzfeld and others, 2024; Liu and others, 2024) in the horizontal direction. However, heavy crevassing of the ice surface characteristic during a surge can complicate elevation analysis at high resolution resulting in more coarse analysis in the vicinity of crevasse fields (e.g. 250 m resolution in crevassed areas in Liu and others, 2024).

In this paper, we present an approach that facilitates high-resolution (1–5 m) analysis of the evolution of geophysical processes during a surge, including height changes, crevassing, mass transfer through the glacier system and ice-surface roughness. We utilize data from two years of ICESat-2 observations collected during the mature phase of the Negribreen surge in 2019 and 2020 to demonstrate the measurement capabilities of the ICESat-2 lidar altimeter to resolve the spatiotemporally complex ice-dynamic processes that occur during a surge.

1.2 The study area: Negribreen, Svalbard

The Negribreen Glacier System, located in eastern Svalbard (Fig. 1), surged in late 2016 reaching speeds of 21m d⁻¹ during its peak in July 2017, equivalent to 200 times its normal quiescent velocity (Strozzi and others, 2017; Herzfeld and others, 2021). In response to this rare event, the authors’ Geomathematics, Remote Sensing and Cryospheric Sciences Group at the University of Colorado, Boulder, conducted three airborne survey campaigns of the glacier system in the summers of 2017, 2018 and 2019 (Herzfeld and others, 2022), whose data supplement the analysis in this paper (see Section 2.2).

Figure 1. ICESat-2 survey lines over the Negribreen Glacier System. The survey lines for each of ICESat-2's three beam-pairs are color coded by their reference ground track (RGT) while the Negribreen Glacier System borders are given by the black line. The thick lines correspond to the part of the track that is analyzed in this paper which is mostly equivalent to the boundaries of Negribreen Glacier. Left/Right (L/R) beam-pairs are separated by ~90 m on-ice which is within line thickness over Negribreen in this figure. The insert in the upper right gives the location of the Negribreen Glacier System (red box) within the Svalbard archipelago. Background image from Landsat-8 acquired 5 August 2019.

The Negribreen Glacier System consists of Negribreen, where the majority of the surge activity occurs, Rembebreen, a southern tributary glacier in the upper glacier system, and two main tributary glaciers flowing in from the north: Akademikarbreen that feeds Negribreen in the upper glacier, and Ordonnansbreen further downglacier (Fig. 1). The Negribreen Glacier System receives large amounts of inflowing ice from the Filchnerfonna accumulation zone above the Negribreen Glacier System to the west. The divide between Filchnerfonna and the Negribreen Glacier System as we have defined it in this paper is somewhat arbitrary as the two ice masses are connected by a series of glaciers and ice falls (e.g. Transparentbreen and Filchnerfallet) and are dynamically connected, which we show in this paper. However, the vast majority of the surge activity, particularly in 2019 and 2020, occurs below the Filchnerfonna and as such we define the bounds of the Negribreen Glacier System using the black line in Figure 1, which excludes Filchnerfonna.

Negribreen consists of polythermal ice and is marine-terminating and thus the mechanisms leading to surge behavior, along with surge evolution as a whole, differ from those of a surge in a temperate and/or land-terminating glacier such as Bering Glacier, Alaska (Dowdeswell and others, 1984; Murray and others, 2003; Trantow, 2020), though some unifying surge theories have been postulated (e.g. Benn and others, 2019). Similar to other observed surges in tidewater glaciers in Svalbard (Strozzi and others, 2017; Nuth and others, 2019), Negribreen began accelerating near the terminus as signified by a collapse in the frontal area (location indicated in Fig. 2a), and surge effects, such as crevassing and increased velocities, propagated upglacier affecting other parts of the greater Negribreen Glacier System (Haga and others, 2020; Herzfeld and others, 2021). High velocities and enhanced calving occurring during the surge led to extensive mass loss in the system, which has contributed to sea-level rise.

Figure 2. Landsat-8 images of the Negribreen Glacier System 2018–2021. Landsat-8 RGB imagery (30 m resolution) acquired (a) 30 September 2018 with the upglacier extent of the collapsed surface that occurred during surge initiation circled in black, (b) 20 August 2019 with black lines indicating the lower part of the Negribreen–Akademikarbreen Medial Moraine (NAMM) which has ‘folded’ due to the surge, (c) 31 July 2020 and (d) 8 August 2021 with Ordonnansbreen's tooth indicated by the black arrow.

The mass loss and disintegration of the glacier system is immediately apparent by examining the time series of Landsat-8 imagery in Figure 2 from 2018 to 2021 (RGB imagery at 30 m resolution, Roy and others, 2014). The Negribreen Glacier System terminus has seen massive deformation with the ‘tooth’ of Ordonnansbreen becoming completely detached from the main system by 2021 (see Fig. 2d). Detailed altimeter documentation of the Ordonnansbreen tooth detachment is provided by the ICESat-2/DDA-ice time series in Section 4.3.3. Additional Landsat-8 imagery of lower Negribreen in 2015–2017 documenting surge initiation is provided in the Supplementary Material (negri.change.suppl.pdf, Fig. S29).

1.3 ICESat-2 and the DDA-ice

NASA's ICESat-2 satellite, launched on 15 September 2018, provides high-resolution height measurements of ice sheets and glaciers via its state-of-the-art micro-pulse photon-counting lidar technology (Markus and others, 2017; Neumann and others, 2019). ICESat-2's payload, the Advanced Topographic Laser Altimeter System (ATLAS), provides height observations at a nominal 0.7 m along-track resolution (under clear-sky atmospheric conditions), a significant upgrade to its predecessor, ICESat's Geoscience Laser Altimeter System (GLAS), which provided height estimates every 173 m along-track from 2003 to 2009 (also constrained by atmospheric conditions) (Zwally and others, 2002).

ATLAS operates at a wavelength of 532 nm (green light) and consists of 3 pairs of strong and weak beams with each pair separated by 3.3 km in the across-track direction (Neumann and others, 2019; Martino and others, 2019). Each strong beam is separated from its weak counterpart by 90 m across-track and 2.5 km along-track. The ratio of transmit energies between strong and weak is approximately 4:1. The transmitter array is rotated approximately every 6 months so that either the weak beams lead (flying ‘forwards’, spacecraft orientation = 1), or the strong beams lead (flying ‘backwards’, spacecraft orientation = 0)(see Fig. 3 and Table 1 of Herzfeld and others, 2021).

ICESat-2 flies at ~500 km altitude and orbits in a 92° inclination, and thus provides measurements up to $\pm$88° latitude. The satellite orbits in an exact 91 d repeat cycle with 1387 reference ground tracks (RGTs) (Neumann and others, 2019, 2020b). All-in-all, ICESat-2 provides six measurement opportunities across a 6.6 km swath every 0.7 m along its flight path, repeating every 91 d.

Several factors limit the ideal spatiotemporal measurement resolution of ICESat-2 in the cryosphere, the most significant being cloud cover which severely attenuates green light. Other factors contributing to weak or non-existent return signals over ice include diffuse scattering, the presence of water or dirt and pre-scheduled satellite maneuvers that either pause active measurement for maintenance or calibration purposes, or point the instrument off its normal track in order to better capture a particular phenomenon (Luthcke and others, 2021; Magruder and others, 2021). Signals are further complicated by noise emanating from the solar background or from the instrument itself (Martino and others, 2019). Furthermore, NASA's official ice-surface height product ATL06 (Smith and others, 2020) only resolves surface heights with 40 m postings at 20 m spacing and does not take full advantage of the high-resolution capabilities of ICESat-2 (Herzfeld and others, 2021).

Geolocation of repeat passes may vary up to 25 m across-track in 2019–2020, as demonstrated for each of the 6 beams in the Supplementary Materials (Figs S14–16), due to pointing errors. Moreover, total position pointing control error was over 1800 m at the start of the mission until a correction was uploaded to ICESat-2 on 28 March 2019 (Luthcke and others, 2021). As a result, height estimates in this analysis from January 2019 through March 2019 are significantly offset from the other repeat passes and are therefore not used in the time series analysis of Section 4.3.

The Density-Dimension Algorithm for ice surfaces (DDA-ice), however, fully exploits ICESat-2's measurement capabilities by identifying height-signals at the sensor resolution of 70 cm (Herzfeld and others, 2017, 2021). At this resolution, important geophysical processes can be measured such as crevassing, calving and rifting (Herzfeld and others, 2021). In the present study, we employ the DDA-ice in order to maximize ICESat-2's ability to capture geophysical processes in a rapidly changing glacier system that is actively surging.

2. Data

2.1 ICESat-2 ATL03 data

The DDA-ice algorithm takes as input the ATLAS/ICESat-2 L2A Global Geolocated Photon Data (ATL03) (Neumann and others, 2020b). We use release 4 (revision 1) of these data in this analysis (Neumann and others, 2020a). The ATL03 data give height above the WGS 84 ellipsoid (ITRF2014 reference frame) along with latitude, longitude and time for all the photons downlinked from ICESat-2.

ATL03 data are segmented by granules that each cover 1/14th of a single orbit. Granules are given in HDF5 format and are freely available through NASA or from the National Snow and Ice Data Center (NSIDC). Granules are named using the format ATL03$\_$[yyyymmdd][hhmmss]$\_$[ttttccss]$\_$[vvv_rr].h5, where[yyyymmdd][hhmmss] is the date and time of acquisition associated with the first data point in the granule, [tttt] is the RGT number, [cc] the cycle number, [ss] the segment number and [vvv$\_$rr] the version and revision numbers.

For our analysis, we attain all granules covering the Negribreen Glacier System between 1 January 2019 and 31 December 2020. With the Negribreen Glacier System lying between latitudes 78.5°N 78.8°N, the segment number is equal to either 03 if the satellite is ascending during data collection or equal to 05 if the satellite is descending. The cycle number identifies the number of 91 d cycles that have elapsed since ICESat-2 entered its science orbit. Our analysis here uses data that span eight cycles from cycle 02 to cycle 09. Each cycle is divided into 1387 unique RGTs. There are nine RGTs that provide significant data coverage of the Negribreen Glacier System: 91, 152, 389, 450, 594, 831, 892, 1036 and 1334. Nine RGTs spanning eight cycles equates to 72 granules for this two-year analysis. With all six beams covering part of the glacier system for each RGT, we process a total of 432 measurement passes in this analysis. The specific granules used in this analysis are identified in the Supplementary Material (negri_data_2019_2020.xlsx).

ICESat-2's six beams are labeled as gt1l, gt2l, gt3l, gt1r, gt2r and gt3r. The strong beams can be associated with either the left (l) or right (r) side depending on the orientation of the ICESat-2 observatory, which switches every 6 months or so in order to maximize the solar illumination of the solar panels. At the beginning of 2019, ATLAS was in its ‘backward’ orientation with left ground tracks corresponding to strong beams and right ground tracks corresponding to weak beams. ATLAS switched to its ‘forward’ orientation in September 2019 resulting in strong beam identification by the right ground tracks. A switch back to the ‘backwards’ orientation occurred in June 2020.

2.2 Airborne altimeter and image data

The authors, along with other members of the Geomathematics, Remote Sensing and Cryospheric Sciences research group, collected airborne imagery and altimeter data in the summers of 2017, 2018 and 2019 as part of a campaign to document the surge of Negribreen (Herzfeld and Trantow, 2021; Herzfeld and others, 2022). Each field campaign gives a synoptic view of the entire glacier at a single point in time, providing a tie-in for satellite data, e.g. ICESat-2 tracks, which are spaced out. The analysis in the current paper utilizes the imagery collected during the 2019 campaign, which under-flew several ICESat-2 tracks in upper Negribreen on 13 August 2019. Airborne altimetry data are used to validate the crevasse spacing and crevasse depth, as shown in Herzfeld and others (2021) and Herzfeld and others (2022), and complement our findings presented later in the paper. Here, we give a quick summary of the 2019 Negribreen campaign and the instruments that were employed to document the surge.

For the 2019 campaign, height measurements were attained using a LaserTech Universal Laser System (ULS) instrument that operates at 905 nm. The ULS was used in conjunction with a 1 Hz LORD 3DM-GX5-15 Virtual Reference Unit, i.e. an inertial measurement unit (IMU), along with a 10 Hz Trimble R10 rover GPS in order to attain accurate height estimates. Both the ULS and GPS were attached to the skids of a helicopter which was flown 100–200 m above the ice surface during operation of the laser. With an effective measurement rate of 400 Hz, the ULS provided glacier surface heights for a single nadir profile every 0.06–0.08 m along-track.

Specifically, the 2019 flight campaign under-flew two beam pairs: RGT 594 (gt1l and gt1r) from 5 August 2019, and RGT 450 (gt1l and gt1r) from 18 August 2019 (Herzfeld and others, 2022). These two tracks passed over interesting crevassed areas in upper Negribreen (Fig. 1) and their dates of collection coincided with the field campaign. Crevasse characteristics along these tracks were quantified in Herzfeld and others (2021) and Herzfeld and others (2022) which found a close agreement between crevasse morphology, spacing and depth between ULS data and ICESat-2 data.

Selected photos from the 2019 campaign, taken with a handheld Nikon D5100 Single-Lens Reflex camera, are given in Figure 3. During the campaign in August 2019, much of the lower glacier was covered in low-lying mixed-phase clouds (Shupe and others, 2011; Gierens and others, 2020) (Fig. 3a). Mixed-phase clouds in Svalbard occur in every season within 1 km of the surface, and have a complicated structure consisting of supercooled liquid and ice layers that obviously obscure the glacier ice-surface (Gierens and others, 2020). These ephemeral low-lying clouds restricted airborne surveys to the upper Negribreen Glacier System for several, but not all, flights in 2018 and 2019.

Figure 3. Imagery from the airborne campaign flights over Negribreen in August 2019. (a) Low-lying clouds covering the lower glacier and terminus (photo looking downglacier). (b) Young surge crevasses in upper Negribreen. (c) Large and complex crevasses exceeding 30 m depth in the center-front of the glacier just upglacier of the terminus. (d) Snow-bridged crevasses seen most clearly in the left-foreground with white, fresh snow covers at the top of the open crevasses. (e) Water-filled crevasses. (f) Crevasses near the Negribreen–Akademikarbreen Medial Moraine (NAMM) that are filled with water indicating a disruption in the local englacial drainage system. Fresh crevassing from mature phase surge dynamics affected this area along the northern NAMM in early 2020, shortly after this photo was taken.

2.3 Sentinel-1 SAR data

We used Synthetic Aperture Radar (SAR) imagery from the European Space Agency's (ESA's) Sentinel-1 satellite (Geudtner and others, 2014) to derive velocity estimates on Negribreen. Mean velocity estimates are derived for several characteristic time periods between 2019 and 2020 which are used to supplement the ICESat-2 analysis by providing basic information on changes in the velocity field associated with the surge.

There are two Sentinel-1 satellites labeled A and B, which in tandem provide repeat imagery every 6 d for a given location. ESA freely provides the Sentinel Application Platform (SNAP) software to derive surface velocity estimates for ice sheets and glaciers using offset tracking methods (Veci and others, 2014). Offset tracking methods measure feature motion between two images using patch intensity cross-correlation optimization.

Over Negribreen, Sentinel-1 operates in the Interferometric Wide (IW) swath mode acquiring data with a 250 km swath at 5 m by 20 m spatial resolution. The SNAP toolbox offset tracking method takes as input the Level-1 Ground Range Detected (GRD) product from two Sentinel-1 images. These images are separated by a temporal baseline, which in our analysis is usually equal to 12 d, i.e. the length of the repeat cycle of a single Sentinel-1 satellite. The GRD products are updated with more accurate orbit data and are then coregistered based on geometry given by the ACE-30 Digital Elevation Model (DEM). We then perform offset tracking to derive a 300 m by 300 m resolution velocity product with missing or removed low-confidence estimates interpolated up to 1 km. A velocity estimate is considered low-confidence if its cross-correlation value falls below a threshold of 0.1 as determined by the SNAP offset tracking software (Veci and others, 2014).

Four characteristic velocity maps from 2019 to 2020 are provided in Figure 4 (see negri_data_2019_2020.xlsx for specific Sentinel-1 data products used). It is clear in each map that Negribreen is surging, reflected by elevated speeds (>1 m d⁻¹), while the tributary glaciers remain at typical quiescent speeds well below 1 m d⁻¹. Velocity magnitudes are largest near the front of Negribreen and decrease across its length with elevated speeds still present near the border with Filchnerfonna in the upper glacier.

Figure 4. Negribreen velocity maps from 2019 to 2020. Maps derived from Sentinel-1 SAR data. (a) Mean surface velocities between 4 February 2019 and 16 February 2019 (m d⁻¹). (b) Mean surface velocities between 11 August 2019 and 23 August 2019 (m d⁻¹). This baseline spans the 2019 airborne campaigns in August 2019. (c) Mean surface velocities between 10 July 2020 and 22 July 2020 (m d⁻¹) with maximum speeds exceeding 8 m d⁻¹. (d) Mean surface velocities between 19 December 2020 and 31 December 2020 (m d⁻¹) at the end of the study period in December 2020 which shows typical velocities during the winter months of 2020. Only velocity estimates with a cross-correlation ratio above 0.1 are displayed (Veci and others, 2014).

Figure 4a displays mean velocity magnitudes near the beginning of our study period between 4 February 2019 and 16 February 2019, with maximal velocities reaching 9 m d⁻¹. Mean velocity magnitudes for Negribreen remain at this level throughout the winter months of early 2019 and do not begin to decrease until August of 2019.

Figure 4b provides mean velocity estimates between 11 August 2019 and 23 August 2019, corresponding to the time of the August 2019 field campaign of (Herzfeld and others, 2022). Maximal surge speeds were near 6 m d⁻¹. As Negribreen transitions into its mature surge phase by 2019, the seasonal component of velocity becomes more apparent as reflected by a significant slowdown in August from maximums typically occurring in July.

Figure 4c provides mean velocities between 10 July 2020 and 22 July 2020 when the glacier system is moving its fastest for the year due to an abundance of meltwater that lubricates the glacier base. Maximal velocities in July 2020 reached 10 m d⁻¹, similar to those in July 2019.

Finally, Figure 4d gives mean velocities between 19 December 2020 and 31 December 2020 near the end of our analysis period. This map provides typical velocities during winter time in Negribreen in 2020 with maximal velocities around 4 m d⁻¹.

3. Methods

3.1 DDA-ice: high-resolution surface heights from ICESat-2 data

We applied the DDA-ice algorithm to the raw photon data found in ATL03 to identify the ice-surface signal at sensor resolution, which is then interpolated at 5 m resolution for smooth ice and 1 m resolution for rough ice, to attain standardized ice-surface height estimates utilized in further analysis. The full mathematical description of the DDA-ice is found in Herzfeld and others (2017) and Herzfeld and others (2021).

The central idea of the DDA-ice is to calculate the density-field for the returned photon point cloud using a convolution defined by a density operator (Herzfeld and others, 2017). Density is then used for classification of signal photons versus background photons in the photon cloud. The density operator, applied to each photon, employs a radial basis function that weights neighboring photons with a 2D anisotropic Gaussian kernel. Signal-noise separation utilizes an auto-adaptive threshold function, i.e a function that automatically adapts to highly variable photon density characteristics including apparent surface reflectance, different background characteristics based on the time of day and some instrument-related artifacts. In addition to identification of surface height in signal photons, the DDA-ice provides interpolated surface heights using a function termed ‘ground follower’ for short. The ground follower employs a piece-wise linear function that is weighted by density values associated with signal photons in the segment of the interpolator. Segment length depends on surface roughness, calculated from signal photons.

The DDA-ice uses algorithm-specific parameters to best identify the type of surface under investigation. Three important parameters that control the shape of the weighting kernel include the standard deviation of the Gaussian distribution (sigma, s), the amount of standard deviations used for weighting (cutoff, u) and the anisotropy factor (a) that specifies the ellipticity of the kernel with positive values giving a more horizontally stretched shape. The two main parameters that control thresholding are the quantile (q) and the offset bias (k). The interpolated surface estimate given by the ground follower has a resolution, R, and provides increased resolution around rough surfaces by a factor of r. The increased resolution is triggered by the roughness parameter, S, that specifies a standard deviation limit for the vertical distribution of thresholded signal photons. Finally, the estimated depth of crevasses is controlled by the crevasse-quantile parameter, Q. The specific parameter values we use in this analysis of Negribreen are given in Table 1.

Table 1. DDA-ice parameters for Negribreen runs in this analysis

The same parameters were used for both the strong and weak beams. Ground-follower resolution when ice surface is rough: With R = 5 and r = 5, then the crevasse-follower resolution is 1 m.

3.2 ICESat-2 ATLAS data processing with the DDA-ice

We ran the DDA-ice on all ICESat-2 ATL03 data over Negribreen from January 2019 through 2020 using the parameters from Table 1, which are close to the default parameters given in Herzfeld and others (2021), but are optimized for the current analysis that spans multiple seasons and uses all the available beams. The same parameter values were used for both the strong and the weak beams. We began by processing all the strong beam data for a given cycle, which consists of three beams per pass with nine passes (individual RGTs) per cycle. After each run we manually filtered out all the cloudy data by looking at the DDA-ice results so as to be sure that we retain only ice-surface signals.

For each beam pass we calculated the along-track ‘ice-signal fraction’ as an indicator of the cloudiness of each measurement pass over the Negribreen Glacier System (see column 2 in Table S1 of the Supplementary Material, negri.change.suppl.pdf). The ice-surface fraction takes the length along-track for which a valid ice-surface height was estimated and divides it by the total along-track survey length. The total along-track survey length (in meters) for each beam pass over the Negribreen Glacier System is given in column 3 of Table S1 of the Supplementary Material (negri.change.suppl.pdf).

Next, we processed all the weak-beam data for which the associated strong-beam pair yielded a non-zero ice-signal fraction. Typically, the weak-beam data provided a non-zero ice-surface fraction when the associated strong-beam data resulted in an ice-signal fraction >0.9.

3.3 Surface height change determination

We analyze ice-surface height change across each of the 54 ICESat-2 ground tracks across Negribreen in 2019 and 2020 in order to estimate mass transfer occurring during the surge evolution. At the resolution of the DDA-ice interpolated ground estimate (1 or 5 m), however, high-resolution morphology, especially crevassing, complicates the bulk height change estimate. For example, fresh surge crevasses in upper Negribreen (Fig. 3b), with depths up to 30 m, lead to an underestimation of the mean ice-surface height, especially when ground-follower resolution is refined over the crevasses. Height change analysis is further complicated if the geophysical signal we are tracking corresponds to crevasse deepening or crevasses advection. Compared to older and wider crevasses however, the cross-sectional area of the crevasse voids for fresh crevasses is relatively small compared to the total width of the glacier. Therefore, for height-change estimates in the case of young surge crevasses in upper Negribreen, we estimate surface height of the glacier surface by using the 90th percentile height for every 30 m along-track bin.

Crevasses in lower Negribreen are more complex and have undergone several deformational processes at this point in the surge. While originally they may have resembled the young crevasses in upper Negribreen, large-scale dynamics near the front of the glacier have transformed these mature crevasses to appear markedly different from those in Upper Negribreen as seen in Figure 3c. Here, the cross-sectional area of crevasse voids is significantly larger with respect to the total glacier width across-track. We therefore calculate surface heights for height-change analysis in this region by taking the 50th percentile surface height for every 30 m along-track bin. This gives the average surface height within the 30 m bin and smooths out processes that change individual crevasse characteristics that may complicate surface height change determination.

3.4 Roughness and crevasse characteristics

From ICESat-2/DDA-ice surface heights, vario functions are calculated to derive surface roughness values, which are characteristic of the spatial structure of the ice surface (Herzfeld, 2008; Herzfeld and others, 2021). Every 200 m along-track, we calculate discrete first-order vario function vectors, v^k = [v ₁,  …,  v _j,  …,  v _N], within 400 m windows centered at location x _k, which act on n pairs of height estimates, z(x _i) and z(x _i + h) separated by some lag distance h grouped in bins defined by the vector h, whose bounds are indexed by j:

(1)$$v_j^k = {1\over 2n} \sum_{i = 1}^{n} \big(z( x_i) - z( x_i + h) \big)^2$$

where x _i is the along-track location of a height estimate (z) and h is the separation, or lag, distance between pairs of points with h _j−1 < h ≤ h _j, for j = 1,  ...,  N, with h ₀ = 0. We use N = 10 discrete lag distance bins of length 40 m to characterize the surface at a given point x _k, implying h = [0,  40,  80,  …,  400].

In cases where there are underlying regional trends in the data, such as a glacier's surface slope, it is more useful to use the residual vario function to quantify roughness (Herzfeld, 2008). Using the mean value m at a center point x _k, given by

(2)$$m_j^k = {1\over n} \sum_{i = 1}^{n} \big(z( x_i) - z( x_i + h) \big)$$

the residual vario function V for a reference center point x _k is defined as

(3)$$V_j^k = v_j^k - {1\over 2} m_j^k$$

where lag distance bins h _j−1 < h ≤ h _j, for j = 1,  …,  N are used in both Eqns (2) and (3).

We derive an estimate of surface roughness, ζ _k, at point x _k by taking the maximum of V^k:

(4)$$\zeta( x_k) = max( {\boldsymbol V^k}) $$

The parameter ζ is equivalent to the 1D pond parameter for a residual vario function (respond) introduced in Herzfeld (2008), though in the current analysis we refer to this quantity simply as roughness (of the glacier surface). This parameter is also derived using 2D (residual) vario functions applied to satellite imagery of Bering Glacier, Alaska, during surge in Trantow and Herzfeld (2018).

Additional crevasse and crevasse field characteristics, i.e. crevasse spacing and depth, are calculated from the DDA-ice height outputs for selected ICESat-2 tracks. Crevasse spacing and depth are determined using a simple deterministic algorithm (Herzfeld and others, 2013b, 2021, 2022) that identifies individual crevasses when the so-called ‘jump height’ between crevasse top and crevasse bottom exceeds 2 m.

3.5 Roughness change determination

We calculate mean rates of roughness change in 2019–2020 using the roughness measure, ζ, derived in Eqn (4). While we are interested in the year-to-year dynamical component of roughness change, it is important to consider its significant seasonal component. During the winter months of snow accumulation, approximately October through April, crevasses of moderate and narrow widths can become covered, or bridged, by snow. These snow-covered crevasses can persist throughout the summer in locations of low deformation as seen in Figure 3d. In this case, only the widest crevasses will be detected by the DDA-ice, which finds the primary signal. The DDA-bifurcate-seaice algorithm (Herzfeld and others, 2023) is able to detect secondary signals such as those beneath the water surface, but is not employed in this analysis as it is not yet reliably adapted to land-ice applications. Because many crevasses are bridged by snow in winter, roughness change is derived using data from summer months only when snow-bridges are the least prevalent (roughly May through October). During peak melt season in July and August however, crevasses exist that are filled with water at locations where the surge deformation has destroyed normal englacial drainage paths that route water from the glacier surface to the base (see Fig. 3e) (Kamb and others, 1985; Harrison and Post, 2003). Water-filled crevasses are present throughout the glacier system but affect only a small number of the total crevasses (Herzfeld and others, 2022), and therefore should not significantly affect the overall roughness estimates in this analysis.

4. Results

4.1 Surface height and surface height change

4.1.1 Surface heights and data coverage for each ICESat-2 cycle

Along-track ice-surface height estimates yielded by the DDA-ice ground follower at 1–5 m resolution, as determined by the algorithmic parameters in Table 1, are given for the four cycles in 2019 and the four cycles in 2020 (Fig. 5). Heights in the Negribreen Glacier System range from ~0 m (sea-surface height) at the terminus to over 800 m in the upper tributary glaciers. The most sparse coverage of the glacier system occurs during the late-summer/autumn months (Figs 5c, g), particularly in the lower part of the glacier system, when low-level mixed-phase clouds have their highest occurrence (Shupe and others, 2011; Gierens and others, 2020) (see photograph in Fig. 3a). During all but a couple cycles, Akademikarbreen is more densely covered than Ordonnansbreen, likely due to the local orography, which strongly influences the occurrence and characteristics of mixed-phase clouds (Gierens and others, 2020).

Figure 5. Negribreen ice-surface elevation estimates provided by the DDA-ice for single ICESat-2 cycles (91 d). 2019 estimates (top row): (a) January–March 2019, (b) April–June 2019, (c) July–September 2019, (d) October–December 2019. 2020 estimates (bottom row): (e) January–March 2020, (f) April–June 2020, (g) July–September 2020, (h) October–December 2020.

In general, clouds are the main factor controlling the ice-signal fraction for each beam pass and the survey coverage for each ICESat-2 cycle. Optically thick clouds will fully attenuate ICESat-2 532 nm transmit photons resulting in ice-signal fractions equal to zero, while optically thin clouds only partially attenuate the transmit energy resulting in ice-signal fractions <1 (but possibly >0). Therefore, the coverage maps shown in Figure 5, along with quantified ice-signal fraction estimates given in the Supplementary Material, provide seasonal cloudiness information that are useful for airborne campaigns such as those conducted by Herzfeld and others (2022).

4.1.2 Surface height rate of change (2019–2020)

Figure 6 gives the average surface height rate of change from 2019 to 2020 following the processes outlined in Section 3.3. While this calculation provides mean height changes in meters per year, it is important to note that the majority of the change occurs in the summer months, with very little height change between January and April. This is seen more clearly in the times series plots of Section 4.3.

Figure 6. Rates of change of glacier surface height during the 2019–2020 part of the recent Negribreen surge. Surface height change rate in meters per year.

Figure 6 shows general thinning in the upper glacier above 5.85 × 10⁵ m UTM-East and general thickening in the lower glacier below, implying a mass transfer from the upper glacier to the lower glacier. Surface heights grew between 2019 and 2020 in the front ~1/3 of the glacier at a rate reaching 30 ma⁻¹ near Negribreen's terminus. A chaotic pattern of tightly spaced surface lowering and surface gains exceeding ±30 ma⁻¹ near the center-front of the glacier reflects the generation and advection of massive and complex crevasses (see Fig. 3c).

Surface lowering occurs in the upper ~2/3 of Negribreen, while surface heights remain mostly constant on tributary glaciers such as Akademikarbreen and Rembebreen, which is seen clearly at the crossing of the medial moraines in Figure 6. Aside from a slight lowering near its terminus, Ordonnansbreen also experiences almost no change in surface height from 2019 to 2020 despite the significant height changes that occur just across the medial moraine on the surging Negribreen.

There is, however, a clear surface lowering of around 8  ma⁻¹ on Transparentbreen, which connects Negribreen to the Filchnerfonna accumulation zone above the Negribreen Glacier System. Thus, the surge of Negribreen is affecting additional parts of the glacier system as ice evacuates to the sea, with surge-induced changes expanding beyond the main glacier, across Transparentbreen and into the large catchment area above the glacier (Filchnerfonna). Transparentbreen appears to be the only tributary glacier in the Negribreen Glacier System experiencing significant height-change during the surge expansion from 2019 to 2020.

4.2 Surface roughness and roughness change

Crevassed regions are indicated by high surface roughness (Herzfeld, 2008; Herzfeld and others, 2014; Trantow and Herzfeld, 2018; Herzfeld and others, 2021, 2022). Negribreen's ice-surface roughness, given by ζ in Eqn (4), for the four ICESat-2 cycles of 2019 and 2020 is shown in Figures 7 and 8 respectively. The medial moraines provide a clear boundary between the rough and surging Negribreen and the relatively smooth surfaces of the non-surging tributary glaciers. In general, roughness is larger further downglacier as the surging ice experiences more deformation events moving through the glacier system with the strongest events experienced near the glacier front where surge velocities are largest (Fig. 4).

Figure 7. 2019 Negribreen roughness data per cycle as given by the natural logarithm of the ζ parameter. (a) January–March 2019, (b) April–June 2019, (c) July–September 2019, (d) October–December 2019.

Figure 8. 2020 Negribreen roughness data per cycle as given by the natural logarithm of the ζ parameter. (a) January–March 2020, (b) April–June 2020, (c) July–September 2020, (d) October–December 2020.

Changes in seasonal roughness characteristics, seen by comparing plots within Figures 7 and 8, must be interpreted with caution as seasonal effects, such as the loss of snow cover, are intertwined with dynamically induced roughness changes (see Section 3.5). With this in mind, we can use these maps to identify timing of significant roughness changes. For example, we see persistent large roughness values along the Negribreen–Rembebreen boundary beginning in the summer of 2019 (lnζ > 1, Fig. 7c), indicating the timing of a large surge-induced deformation event in the area. Similarly, we see a gradual increase in surface roughness at the Negribreen–Akademikarbreen boundary with roughness estimates (lnζ) growing from well below 0 in 2019 (Fig. 7), to eventually surpassing 1 by the summer of 2020 (Fig. 8c). The persistent roughness growth in this area reflects the continued surge activity that occurs in upper Negribreen throughout 2019–2020.

Locations of surge activity are more clear in Figure 9 which gives the mean rate of change of dynamically induced roughness over the 2019–2020 period of the Negribreen Glacier System surge following the approach discussed in Section 3.5. Recall that Negribreen began accelerating in 2016 near the glacier front and surge-induced changes have expanded upglacier as the surge phase matures (Haga and others, 2020; Herzfeld and others, 2022). In Figure 9 we see that for the most part, surface roughness has decreased in the lower glacier reflecting reduced surge activity at that location. There are locations in the lower glacier however, that have seen increased roughness indicating the continued occurrence of strong, but isolated, surge deformation events.

Figure 9. Rates of change of roughness during the 2019–2020 part of the recent Negribreen surge. Roughness (ζ) change rate in Δζ per year.

In contrast, large portions of the upper glacier experienced increased roughness over the study time interval (e.g. red regions in Fig. 9), which reflects expansion of surge effects upglacier in Negribreen. Large increases in surface roughness also occurred on the ice falls between the Filchnerfonna and Negribreen, which illustrates further expansion of surge effects beyond Negribreen and into the accumulation zone above the glacier.

4.3 Time series of ICESat-2 surface height profiles

In this section, we analyze time series of selected ICESat-2 profiles that allow derivation of glaciological changes associated with the surge of the Negribreen Glacier System in 2019–2020. We present a subset ICESat-2 profiles, identified by their associated RGT and beam, based on the glaciological insight that their result provided. The full collection of ICESat-2 time series for each of the 54 ICESat-2 profiles can be found in the Supplementary Material (negri.change.suppl.pdf).

4.3.1 Detection of new surge crevasses and changes in existing crevasse fields

ICESat-2 data analyzed with the DDA-ice facilitate identification of crevasses (Herzfeld and others, 2017, 2021, 2022). Here we utilize this crevasse detection capability to analyze and map the progression of the surge through the Negribreen Glacier System in 2019–2020.

All crevasse fields observed and analyzed here evolved morphologically through expansion or contraction, advection downglacier, formation of a snow-bridge and/or filling with water. An example of crevasse expansion is given by the occurrence of new crevasse fields formed in upper Negribreen along RGT 594 gt1l as seen in Figure 10a between 8.727 and 8.728 × 10⁵ UTM-North where a smooth, uncrevassed ice-surface is reported in early 2019 (orange and yellow lines), and at the same location new crevasses have opened by November 2019 (green line). This particular finding indicates that additional crevasse fields formed in the southern part of upper Negribreen in late 2019, south of an existing larger crevasse field closer to the boundary with Akademikarbreen. Average crevasse spacing and maximal depth of crevasses along RGT 594 gt1l near Akademikarbreen were estimated using airborne laser altimeter measurements during the field campaigns in 2018 and 2019, which matched estimates derived from ICESat-2 DDA-analyzed data (Herzfeld and others, 2021, 2022). Another example of the detection of crevasse expansion is found along RGT 594 gt2l in Figure 10b near 8.726 × 10⁵ m UTM-North where large crevasses formed along Negribreen's southern margin between February and April 2020.

Figure 10. DDA-ice results for RGT 594 over a surge affected area in upper Negribreen in 2019–2020. (a) RGT 594 gt1l, (b) RGT 594 gt2l and (c) RGT 594 gt3l.

Further evidence of the expansion of surge effects in upper Negribreen is given by the time series of RGT 450 in Figures 11a,b. Large crevasses reaching 10 m depths near the southern margin (left-side of figure) formed in early 2020 similar to those detected in the nearby RGT 594 data (Fig. 10). Note that beam-pair 1 of RGT 450 is also surveyed and analyzed in Herzfeld and others (2022), which shows consistent detection of ~10 m deep crevasses in both the ICESat-2/DDA-ice and the airborne laser altimeter data in 2019. In addition, Figure 11c, which shows RGT 152 gt3r, shows crevasse expansion across both the northern (~8.732 × 10⁵ UTM-North) and southern margin (~8.727 × 10⁵ UTM-North) between early 2019 (orange line) and 2020 (green, blue, brown lines). RGT 152 gt3r also presents an example of isolated crevassing and crevasse enlargement in the center of upper Negribreen near 8.729 × 10⁵ m UTM-North.

Figure 11. Crevasse evolution in upper Negribreen as given by ICESat-2's RGT 450 and 152, 2019–2020. (a) RGT 450 gt1l, (b) RGT 450 gt1r and (c) RGT 152 gt3r.

The apparent disappearance of crevasses from 2019 to 2020, such as those in Figure 10c near 8.729 × 10⁵ m UTM-North, may result from crevasse contraction, advection downglacier without upglacier replacement of similar crevasses, snow-bridge formation or filling with water. Imagery may be used to determine the exact process involved. The airborne photograph in Figure 3a, taken near 8.729 × 10⁵ m UTM-North, shows that the crevasses under consideration, i.e. those found in and to the left of the medial moraine in the image, have become bridged by snow.

Returning to RGT 594 gt1l, we proceed to investigate changes in the width and depth of crevasses in the large crevasse field seen in Figure 10a between 8.730 and 8.732 × 10⁵ UTM-North. Mean surface height change (surface lowering) along this 2 km segment was −11.66 m between August 2019 and August 2020. Mean crevasse spacing increased from 55.45 to 57.6 m, mean depth decreased from 10.71 to 9.82 m, and maximum depth decreased from 16.01 to 15.82 m. These changes in the crevasse characteristics are typical of changes in crevasse fields that formed during earlier years of the surge in 2017 or 2018, as recorded during the field campaigns. This result is consistent with the plot of Figure 9 that reveals this region to be the only one in upper Negribreen that saw a decrease in surface roughness from 2019 to 2020. Over time, crevasses tend to widen a little, become shallower, and the crevasse edges are more rounded, as erosion progresses (Herzfeld and others, 2013b).

In general, we find the local maximum of crevasse depth increases from upglacier to downglacier regions of Negribreen. From ~16 m depth along RGT 594 gt1l in upper Negribreen, maximal depths increase to ~22 m near the mid-glacier along RGT gt3l (Fig. 10c) and exceed 30 m in the lower glacier near the terminus.

4.3.2 Crevasse expansion along and across the shear margin

The Negribreen–Akademikarbreen Medial Moraine (NAMM) is easily identified in airborne and satellite imagery as the dark dividing line between the surging ice of Negribreen and the non-surging ice to the north. From its formation point at the Negribreen–Akademikarbreen junction, the NAMM advects downglacier, past the Lykkenhøgda hills at mid-glacier and along the Negribreen–Ordonnansbreen border in the lower glacier, leaving an obvious stripe through Negribreen near its northern margin. The NAMM provides an example of the folded moraine that can be used to identify a glacier as a surge-type glacier (Post, 1972; Lefauconnier and Hagen, 1991). The fold results from a shift in the dynamic equilibrium between a surge-type glacier (here, Negribreen) and a neighboring, non-surging glacier (Ordonnansbreen). The evolution of the lower NAMM, along with its folds, is seen clearly in the Landsat-8 imagery in Figure 2.

Throughout Negribreen's longitudinal extent, large shear stresses exist around the NAMM due to its boundary between the surging ice of Negribreen and the non-surging ice of the northern tributary glaciers. In 2017 and 2018, Herzfeld and others (2022) documented the disintegration of the lower NAMM, i.e. the Negribreen–Ordonnansbreen shear margin, along which the terminus retreated via processes of rifting and calving at a pace faster than observed elsewhere in the glacier (visualized as a ‘retreating bay’ at the terminus). The formation of bays and melange areas occur at former areas of so-called chaotic crevasse types (Herzfeld and others, 2013a), which are often found at locations with large shear stresses.

In the current analysis of surge progression, we detect the development of additional large crevasses in northern Negribreen in the upper-mid glacier around the NAMM (i.e. upglacier of the largest ‘retreating bay’ observed in 2017) in 2020. Figure 12 displays a transverse expansion of crevassing (outward toward the NAMM), recorded in RGT 892 gt1l and gt1r. In August 2019 (yellow lines) deep crevasses (~10 m) are measured across the majority of the glacier with a northern extent around 8.729 × 10⁶ m UTM-North. By August 2020 (purple lines), 1 year later, these large 10 m deep crevasses expand transversely (northward) ~5 km to the extent of 8.7295 × 10⁶ m UTM-North. These new crevasses appear to be wider than the older ones to the south and are an indication of the intensification of the shear margin.

Figure 12. DDA-ice results near the Negribreen–Akademikarbreen Medial Moraine (NAMM), 2019–2020. (a) RGT 892 gt1l and (b) RGT 892 gt1r.

In August 2019, just before this noted crevasse expansion, large amounts of water were observed in the crevasses near the northern NAMM as seen in Figure 3f. The presence of clear blue meltwater in crevasses is indicative of a local disruption in the englacial drainage system, which occurs in an actively surging region (Kamb and others, 1985; Trantow and Herzfeld, 2024).

In Figure 13, gt2l and gt2r of RGT 1334 show the same transverse (northern) expansion of crevassing across the NAMM at the mid-glacier. The large 10 m deep crevasses expand toward the NAMM beginning at 8.7288 × 10⁶ m UTM-North in September 2019 (yellow lines), up to 8.7289 × 10⁶ m UTM-North by June 2020 (blue lines) and to 8.729 × 10⁶ m UTM-North by December 2020 (brown line, visible only in gt2l due to partial cloud cover). These changes are reflected in Figure 9 where increased roughness is observed from 2019 to 2020 at locations near the NAMM in upper Negribreen. Additional evidence of crevasse expansion along and across the NAMM is given in Figure 11c for RGT 152 gt3r which surveys Negribreen just upglacier of RGT 1334 gt1l/r.

Figure 13. DDA-ice results near the Negribreen–Akademikarbreen Medial Moraine (NAMM), 2019–2020. (a) RGT 1334 gt2l and (b) RGT 1334 gt2r.

Our ICESat-2 results of fresh crevassing in mid and upper Negribreen around the NAMM indicate an expansion of shear margin disintegration upglacier from a retreating bay at the terminus. Taken together, the NAMM disintegration and the retreating bays at the terminus illustrate the evolution of a shear margin with a very strong velocity gradient along a folded moraine (seen most starkly along the lower NAMM in the velocity maps in Fig. 4). This result demonstrates how time series of ICESat-2 data, analyzed with the DDA-ice, can be employed to derive quantitative information about complex surge processes, here, transverse deformation across a medial moraine forming a shear margin with increasingly large stress gradients as the surge evolves. To our knowledge, this is a new capability of satellite-based observation.

4.3.3 Disintegration of the Ordonnansbreen tooth

The Ordonnansbreen tooth, a small area of ice indicated by a black arrow in Figure 2d, has been detaching from Ordonnansbreen's north terminus since the surge began in 2016. As noted earlier, Ordonnansbreen has been moving at quiescent speeds throughout our total observation time (2016–2021), as seen in the velocity maps of Figure 4.

The detaching tooth, along with the greater disintegration of the ice along the nearby medial moraine (lower NAMM), is visualized in the time series of Landsat-8 RGB imagery from 2018 to 2021 in Figure 2. In 2018 (Fig. 2a) the tooth is attached to the northern part of the Ordonnansbreen terminus. The tooth's connection with the terminus shrinks by August 2019 (Fig. 2b) and becomes totally disconnected by July 2020 (Fig. 2c). Ordonnansbreen's terminus rapidly retreats between July 2020 and August 2021 (Fig. 2d), greatly increasing the distance between the shore-fast tooth and the calving front. While the Landsat-8 imagery helps visualize this process in two spatial dimensions, only with the ICESat-2 data do we get surface height and rifting-depth information to elucidate the third dimension of this rapidly changing glacial feature.

We see rifting in the vicinity of Ordonnansbreen's tooth in the time series of RGT 91 gt3l and gt3r (Fig. 14), which both cross the terminus of Ordonnansbreen. The orange line in Figure 14, representing the ice-surface in April 2019, shows continuous ice extending from the shore line at 8.73 × 10⁶ m UTM-North to the calving front near 8.7274 × 10⁶ m UTM-North. By January 2020, a ~10 m deep rift developed near 8.7299 × 10⁶ m UTM-North that penetrated all the way to the ocean surface, disconnecting the glacier from the shoreline (see cyan-colored line in Fig. 14). The rift geometry at this location remained relatively fixed between January 2020 (cyan) and July 2020 (purple, gt3l only) indicating little surge activity during this time. By September 2020 (brown lines), the rift between Ordonnansbreen and the shoreline increased significantly from ~400 m to over 1.5 km in the direction of the RGT 91 survey line. The September 2020 result in RGT 91 gt3r (Fig. 14b) in particular, indicate the presence of large icebergs and an even greater separation distance nearing 5 km between the shoreline and the main glacier. These results indicate that the process of tooth detachment, along with additional nearby disintegration, occurs most rapidly during peak glacier velocities in late summer implying that the initial rift likely developed sometime around July or August 2019.

Figure 14. DDA-ice results over Ordonnansbreen's terminus and the Ordonnansbreen tooth, 2019–2020. (a) RGT 91 gt3l shows the disintegration of Ordonnansbreen's terminus at its northern edge, i.e. at the tooth, above 8.729e6 UTM-North. (b) RGT 91 gt3r also displays the signal of a disintegrating northern terminus of Ordonnansbreen. Both time series also show significant surface lowering across the terminus width. Note that the Landsat-8 image in the insets used for visualization of track locations is from 5 August 2019.

Looking at the entirety of the RGT 91 gt3l/r survey lines, we observe an ice-surface lowering of 4–5 m across the Ordonnansbreen terminus during the 15 months between April 2019 and July 2020, i.e. ~20$\%$ of its total thickness, indicating rapid height and mass loss. This ice-surface continued to lower by ~1 m across the entire width during the 3 months between July and September 2020. Therefore, while the majority of Ordonnansbreen remains in a quiescent state, the surging Negribreen has a significant effect on the evolution of its neighboring tributary glacier through dramatic surface lowering, rifting and calving at the terminus. The surge effects here, along with the major deformation along its medial moraine with Negribreen (lower NAMM), leave the future state of Ordonnansbreen uncertain as the Negribreen Glacier System surge continues to progress.

4.3.4 Changes in the ice falls between the Filchnerfonna and upper Negribreen

The ice falls between Filchnerfonna and upper Negribreen have been areas of pervasive change during the 2019–2020 part of the surge in the Negribreen Glacier System, as observed during the 2019 field campaign (Herzfeld and others, 2022). While ice falls are characterized by heavy crevassing due to steep topography, the imprint of the surge manifests through fresh crevasse openings along with the widening and deepening of existing crevasses. Figures 15 and 16 show distinct signs of activity in the ice falls, indicative of expansion of the surge-affected area beyond the Negribreen Glacier System. Increased crevassing is detected in the southern ice falls (the largest of which is Filchnerfallet), particularly in 2020, but is not as prevalent along Transparentbreen to the north (Fig. 15). There is less crevassing on Transparentbreen because it is less steep than the other inflowing glaciers from Filchnerfonna, however, there is a pronounced surface lowering across its entire width (Fig. 16). The rate of surface-height lowering here is increasing as seen by the spacing of the later 2020 lines (purple, brown) compared to the earlier 2019 lines (yellow, green), indicating building surge effects on Transparentbreen. In contrast, Filchnerfallet does not experience a clear signal of surface lowering until late 2020, which indicates delayed activation of ice flowing in from southern Filchnerfonna with respect to northern inflowing ice through Transparentbreen.

Figure 15. DDA-ice results near the inflow from Filchnerfonna, 2019–2020. (a) RGT 450 gt2l surveying both the southern ice falls (Filchnerfallet) and Transparentbreen in the north and (b) RGT 450 gt2r surveying the ice falls (Filchnerfallet) and Transparentbreen.

Figure 16. DDA-ice results near the inflow from Filchnerfonna, 2019–2020. (a) RGT 450 gt3l surveying only Transparentbreen and (b) RGT 450 gt3r also surveying only Transparentbreen.

In the context of surge expansion beyond the Negribreen Glacier System, our results motivate the question of whether the surge in the Negribreen Glacier System may induce a disintegration of the surrounding glacial area, potentially destabilizing the Filchnerfonna. A less dramatic interpretation is that significant surface lowering in the Negribreen Glacier System leads to draw-down of ice flowing through the ice falls. Either way, the surge in the Negribreen Glacier System is affecting adjacent ice areas at the border of its primary accumulation zone.

4.3.5 Across-track variability between strong and weak beam pairs

To our knowledge, the DDA-ice is the only algorithm that can provide near-sensor resolution surface heights from both weak beam and strong beam data (Herzfeld and others, 2021), thereby facilitating investigation of the types of spatial changes that persist across the 90 m separation of the two beams in a pair. In the example of RGT 450 gt1l and gt1r (Figs 11a,b), the older crevasse field between 8.7310×10⁶ m UTM-North and 8.725×10⁶ m UTM-North appears similar in both beams. However, the freshly opened crevasse field on the medial moraine in upper Negribreen that separates fast and slow flowing ice (8.728 − 8.729 × 10⁶ m UTM-North) is characteristically different in the two beams. RGT 450 gt1r, surveying 90 m further upglacier, captures a series of crevasses with similar depth whereas RGT 450 gt1l captures a particularly large and wide crevasse (or rift) that lies in the moraine that separates two kinematic provinces. In particular, the acceleration and drawdown of ice flowing in through Transparentbreen to the north of the medial moraine lead to a discontinuity with the slower flowing ice from Filchnerfallet to the south (as seen in Fig. 15). The rift in the medial moraine, detected only in RGT 450 gt1l, shifts southward as the surge progresses, i.e. outward in the transverse direction, which indicates the timing of dynamic and morphological changes that extend the surge-affected area across a larger part of the glacier.

Figure 12 gives a second example of the across-track variability between beam pairs. The crevasse field around 8.729 × 10⁶ m UTM-North is characteristically different in August 2020 (purple line) and November 2020 (brown line) between the two beams, with the strong beam (gt1l) detecting more and larger crevasses than the weak beam (gt1r) that lies 90 m upglacier. This result points to the large spatial variability in crevasse characteristics along the NAMM which divides fast and slow moving ice.

In summary, analysis of across-track variability between ICESat-2 beam pairs provides spatial information at the 90 m length scale, which is on the order of individual crevasse fields (Herzfeld and others, 2021). Therefore, this information can be complementary to image classification of crevasse provinces using machine learning approaches, e.g. Herzfeld and others (2024). In addition, this type of analysis can also inform design of beam geometries for future high-resolution lidar missions.

4.3.6 General mass transfer observations

The surface height changes of 2019–2020 given by Figure 6 reflect significant mass transfer from the reservoir area in upper and mid Negribreen downglacier to the receiving area in the lower glacier within 5 km of the terminus. The ice-mass transferred to the receiving area is eventually transported to the Arctic Ocean via heavy calving during the surge. These mass transfer observations are further detailed in the analysis of surface height time series of the preceding sections.

In total, the ICESat-2 results give a clear indication that the mass transfer and mass loss are surge-induced rather than climatically induced. As seen most clearly in Figures 10a and 11b, surface height profiles overlie each other at the beginning and end of each track (i.e. on the non-surging tributary glaciers on the sides of Negribreen) but significant height change is apparent along the interior of the profile. If surface lowering were climatically caused (mass loss through melting), then similar height changes between 2019 and 2020 would be observed across the entire profile.

5. Conclusions

In this paper, we have derived information on geophysical processes that occurred during the surge of Negribreen Glacier System based on analysis of two years of ICESat-2 ATLAS data from 2019 and 2020, processed with the DDA-ice. ICESat-2 data, analyzed with the DDA-ice, provide a unique and novel capability to obtain geophysical information on high-resolution height changes during a glacier surge from spaceborne altimeter observations. The Negribreen Glacier System provides an ideal study region, as the active surge phase overlaps with the observation phase of the ICESat-2 mission. Negribreen's surge started in 2016, with peak acceleration in 2017 and elevated ice-surface velocities continuing through 2024, while ICESat-2 was launched on 15 September 2018. Here we have analyzed all ICESat-2 data collected between January 2019 and December 2020 over Negribreen to demonstrate our approach. The analysis is supplemented by velocity maps from Sentinel-1 SAR imagery and airborne data from our August 2019 campaign to Negribreen.

Geophysical information on the evolution of the surge is derived from the ICESat-2/DDA-ice high-resolution data, including crevassing, height changes, mass transfer toward the terminus and roughness changes indicative of evolving crevasse fields. Height-change rates in 2019–2020 range from −30 ma⁻¹ in the reservoir areas of upper Negribreen, to +30 ma⁻¹ in the receiving area in the lower glacier near the terminus. Roughness change maps indicate an expansion of the surge in upper Negribreen, particularly near the shear margins, while surge activity in the lower glacier lessened from 2019 to 2020.

Time series analysis of ICESat-2 profiles, analyzed with the DDA-ice, indicates formation of new crevasse fields and expansion of existing crevasse fields, as ice evacuates the system and the surge progresses to affect larger areas of the Negribreen Glacier System. The increased surge activity from 2019 to 2020 in the upper glacier is especially seen on the inflowing glaciers from the Filchnerfonna accumulation zone, and along the NAMM, which divides the surging ice of Negribreen and the non-surging ice of the northern tributary glaciers. This later observation, together with the outward crevasse expansion observed along the southern moraine that separates fast and slow ice, implies an increased zone of actively surging ice across the glacier width as the surge progresses into its mature state.

The fresh surge crevassing along and across the NAMM in the mid and upper glacier in 2019–2020 indicates disintegration along the shear margin and reflects a continuation of the deformational process that resulted in a ‘retreating bay’, an area of open water filled with melange that formed as a result of a strong force gradient along the shear margin between the surging ice of Negribreen and non-surging ice of neighboring Ordonnansbreen, near the terminus in 2017. Occurrence of water in crevasses that reaches almost to the ice surface height was observed near the upper NAMM in August 2019, just before the new crevasses were formed. Furthermore, our analysis of surge progression in 2019–2020 provided detailed height information on the disintegration process of a segment of lower Ordonnansbreen (‘Ordonnansbreen tooth’) from the main glacier along a rift that had formed during the surge.

Large-scale mass transfer, extensive roughness changes and the striking disintegration along the NAMM, and of Ordonnansbreen's tooth and the lower glacier as a whole, clearly illustrate the rapid surface change and mass loss a glacier can experience during a surge. These results demonstrate the novel capability of ICESat-2, analyzed with the DDA-ice, to provide high-resolution height data capable of documenting these complex surge processes from space.