Seasonal ice dynamics in the lower ablation zone of Dagongba Glacier, southeastern Tibetan Plateau, from multitemporal UAV images

Abstract

Most glaciers on the Tibetan Plateau have experienced continuous mass losses in response to global warming. However, the seasonal dynamics of glaciers on the southeastern Tibetan Plateau have rarely been reported in terms of glacier surface elevation and velocity. This paper presents a first attempt to explore the seasonal dynamics of the debris-covered Dagongba Glacier within the southeastern Tibetan Plateau. We use the multitemporal unoccupied aerial vehicle images collected over the lower ablation zone on 8 June and 17 October 2018, and 13 May 2019, and then perform an analysis concerning climatic fluctuations. The results reveal that the mean surface elevation decrease of the Dagongba Glacier during the warm season ($2.81\pm 0.44$ m) was remarkably higher than the cold season ($0.72\pm 0.45$ m). Particularly notable glacier surface elevation changes were found around supraglacial lakes and ice cliffs where ice ablation rates were $\sim$3 times higher than the average. In addition, a larger longitudinal decline of glacier surface velocity was observed in the warm season than that in the cold season. In terms of further comparative analysis, the Dagongba Glacier experienced a decrease in surface velocity between 1982–83 and 2018–19, with a decrease in the warm season possibly twice as large as that in the cold season.

Introduction

Glaciers are important components of global water resources and key indicators of climate change. From 1951 to 2009, the temperature recorded in China increased by 0.24$^\circ$C decade$^{-1}$ with the Tibetan Plateau increasing by $0.3$–$0.4^\circ$C decade$^{-1}$ (Qin and others, 2009; Chen and others, 2015; Liu and others, 2015). As a result of global warming, the glaciers in High Mountain Asia have shown considerable mass losses, especially on the southeastern Tibetan Plateau at $-4.0\pm 1.5$ Gt a$^{-1}$ ($-0.62\pm 0.23$ m w.e. a$^{-1}$) (Brun and others, 2017). The accelerated mass losses of glaciers are generally accompanied by an increase in snowline altitude, thinning of glacier thickness and decrease in glacier velocity (Ruiz and others, 2008; Moore and others, 2009; Neckel and others, 2017; Sakai and Fujita, 2017; Dehecq and others, 2019). In addition, considerable mass losses of glaciers can induce serious natural hazards (Richardson and Reynolds, 2000; Kääb and others, 2005), such as landslides, debris flows (Stoffel and Beniston, 2006; Chiarle and others, 2007) and glacial lake outburst floods (Cenderelli and Wohl, 2003; Ng and others, 2007). Understanding glacier dynamics, especially mass balance and surface velocity, is therefore important for assessing the impact of climatic fluctuations on glacier changes (Cuffey and Paterson, 2010; Paul and others, 2015).

Monsoonal temperate glaciers are characterised by higher accumulation/ablation and mass-balance amplitude, and thus are more sensitive to warming than cold-based glaciers (Brun and others, 2017; Sakai and Fujita, 2017; Zhu and others, 2018). As one of the representative areas of monsoon temperate glaciers in western China, most glaciers on Mount Gongga, southeastern Tibetan Plateau, have been in continuous negative mass balance in recent decades (Su and Shi, 2000; Cao and others, 2019). Quantitative evaluation of their dynamic evolution, particularly seasonal glacier dynamics, is therefore essential for a comprehensive understanding the status of glacier ice in High Mountain Asia and the associated drivers and mechanisms. The Qinghai–Xizang (Tibetan) Plateau Scientific Expedition Team of the Chinese Academy of Sciences conducted systematic field observations of glacier thickness, temperature and movement using radar thickness gauges, thermistor thermometers and steam drills on glaciers in this region from 1982 to 1983 (Li and Su, 1996). Pan and others (2012) demonstrated that the glaciers around Mount Gongga had shrunk by 11.3% ($\sim$2.6% decade$^{-1}$) in area since 1966, and the rate of area loss on the western slope (14.6%) was slightly faster than that on the eastern slope (9.8%). In addition, the lower portions of many temperate glaciers in this region are covered by supraglacial debris (Li and Su, 1996; Xie and Liu, 2010; Neckel and others, 2017). Theoretically, the thick moraine covering the glacier surface can slow somewhat the effects of ablation, but a thin or patchy debris layer has the potential to increase ice melt rates (Østrem, 1959; Li and Su, 1996; Shi and Liu, 2000; Xie and Liu, 2010). It is noteworthy that the impact of debris on glacial ablation and ice dynamics is still not fully understood. Recent studies have provided evidence that the presence of supraglacial lakes and ice cliffs on debris-covered glaciers is associated with rapid and spatially highly variable melting (Sakai and others, 2000; Benn and others, 2012; Immerzeel and others, 2014; Miles and others, 2016; Brun and others, 2016; Salerno and others, 2017; Watson and others, 2017; Brun and others, 2018), which complicates the spatial distribution and seasonal patterns of glacier surface melting.

Conventional field methods comprise ablation stake readings and position measurements, which are generally realised at limited numbers of discrete sites on reachable glacier regions. These detailed field surveys on mountain glaciers are particularly time-consuming and challenging. Currently, space-borne remote-sensing images have been widely used for glaciological research (Scherler and others, 2008; Quincey and others, 2009; Dehecq and others, 2015; Neckel and others, 2017; Altena and others, 2019; Dehecq and others, 2019). Low-cost unoccupied aerial vehicle (UAV) surveys have the advantages of flexible acquisition and highly detailed images, which can overcome many of the difficulties associated with satellite remote sensing and produce accurate digital orthophoto maps (DOMs) and digital surface models (DSMs) as well as textured 3-D models. This method has become an important technology for glaciological research (Hugenholtz and others, 2013; Immerzeel and others, 2014; Bhardwaj and others, 2016; Brun and others, 2016; Watson and others, 2017; Wigmore and Mark, 2017; Brun and others, 2018; Benoit and others, 2019).

In this study, we conducted three repeat low-cost UAV surveys (June 2018, October 2018 and May 2019) of the Dagongba Glacier, an extensive debris-covered temperate glacier located on the western slope of the Mount Gongga. The objectives of this study are twofold: (1) to assess the seasonal (warm and cold) dynamics of glacier downwasting/ablation and ice flow, as well as to discuss the evolution of supraglacial lakes and ice cliffs based on repeated high-resolution mapping, and their impact on mass loss and (2) to investigate the long-term slowdown of the Dagongba Glacier during warm and cold seasons and discuss their driving mechanisms.

Study area

Mount Gongga (29.60$^\circ$ N, 101.88$^\circ$ E) is located in the centre of the transition zone from the Sichuan Basin to the southeastern Tibetan Plateau. It has one of the highest concentrations of modern glaciers in the Hengduan Mountains. The Dagongba Glacier is the longest valley glacier on the western slope of Mount Gongga (Fig. 1), with a length of 10.6 km and an area of 19.1 km$^2$ (Su and others, 2002; Zhang and others, 2008, 2019). The snowline of the Dagongba Glacier is 4900–5100 m a.s.l., with the terminus of its ice tongue at $\sim$3900 m a.s.l. and a slope of <$20^\circ$. In this research, the range between 3900 and 4300 m a.s.l. (covering 4.76 km$^2$ of the ablation area) is investigated to focus on the glacier surface elevation and velocity dynamics.

Fig. 1. (a) Location of the Dagongba Glacier on the western slope of Mount Gongga, southeastern Tibetan Plateau (background image: © Google Earth, 5 October 2017). Showing the entire glacier boundary (the second glacier inventory, Guo and others (2015), our study area boundary (yellow shading) and UAV launch site. (b) The screenshot of the textured 3-D model of the glacier tongue on 17 October 2018.

The ablation zone of the Dagongba Glacier is covered with heterogeneous surface debris, with a local thickness up to 1 m (Cao, 1988; Su and others, 1996; Liu and Zhang, 2017). The debris-covered ice tongue contains scattered supraglacial lakes and ice cliffs, resulting in a complex spatial distribution and higher variability of its surface geomorphology (Li and Su, 1996; Liu and others, 2011). The Dagongba Glacier is a monsoonal temperate glacier and is therefore characterised by high accumulation and intense melting. According to field monitoring by ablation stake readings and positions during scientific expeditions between 1982 and 1983, surface velocities of up to 34.8 m a$^{-1}$ and a maximum ablation of 5.9 m were observed in the upper part of the glacier tongue (Li and Su, 1996). Based on the records from the Gonggasi meteorological station ($\sim \!2$ km, 3700 m a.s.l.), the temperature and precipitation both reach their peak values in summer, with 79% of annual precipitation occurring from May to October (annual sum of 1057 mm) with maximum mean daily air temperature of $13.3^\circ$C, during which the surface velocity of the glacier is 2–3 times that in winter (Cao, 1988; Su and others, 1996).

Data and methods

UAV surveys and data processing

In this study, the Dagongba Glacier was surveyed using an UAV three times on 8 June 2018, 17 October 2018 and 13 May 2019 (Table I). According to the regular pattern of glacier ablation (Li and Su, 1996), our observations could be divided into two seasons: warm and cold. The warm season (8 June to 17 October 2018, a total of 130 d) is accompanied by a higher rate of glacier accumulation and melting, while during the cold season (17 October 2018 to 13 May 2019, a total of 208 d), the nourishment and dynamics of the glacier are relatively weakened.

Table 1. UAV photogrammetry parameters used during the three field campaigns

The surveys were conducted with a fully autonomous rotary-wing UAV of type Phantom 4 Advanced from DJI, equipped with a 20-megapixel CMOS sensor and a prime lens with a focal length of 24 mm (equivalent to 35 mm). For the requirements of photogrammetric processing, the UAV flew autonomously along predefined flight routes designed on-site using DJI GS Pro (Ground Station Pro) software, and the three flight strips were uniformly designed with a normal forward overlap and lateral overlap of 85 and 55%, respectively. The UAV is capable of $\sim$25-min flights in field conditions. Thus, the coverage of the full area of interest requires seven or eight separate flights per session during the three separate expeditions. To ensure the accuracy of subsequent data processing, the flight height of each sortie is designed to maintain a uniform altitude during each operation. During the first data acquisition (8 June 2018), the area marked as bad data in Fig. 2a was heavily obscured by fog raising from the valley, and was masked out during subsequent processing to ensure the accuracy of the results.

Fig. 2. (a) DOM of the Dagongba Glacier produced from UAV imagery acquired in October 2018, the extent of glacier area (dashed black line), area with bad data (due to fog partially obscuring the domain) from 8 June 2018, is masked, GCPs for co-registration (blue hollow rhomb), and (b, c, d) images of some GCPs.

ContextCapture software (Casella and others, 2020) is employed to separately achieve the photogrammetry procedure of the images acquired by the UAV on each acquisition dataset. During the processing, first, the original pictures are oriented by bundle adjustment, and then DSMs and DOMs are generated for three acquisition dates. As the only geolocation information recorded by the onboard GPS is used for aerial triangulation, the initial georeferencing should be calibrated to a unified geographic reference.

There are a number of successful applications of ground control points (GCPs) on mountain glaciers to optimise the photogrammetric processing and assess the uncertainty of the final DOMs and DSMs (Gindraux and others, 2017; Rossini and others, 2018; Groos and others, 2019; Van Tricht and others, 2021). However, because of the rough and unstable surface of the debris-covered Dagongba Glacier, it is challenging to achieve the layout and coordinate collection of GCPs in situ. Therefore, all the subsequent photogrammetric products are co-registered with the 17 October 2018 acquisition to improve the geographic consistency of the three phases of products (Benoit and others, 2019). For this purpose, 13 boulders in the glacier-free area were identified from the bundle adjustment of 17 October, and used as the virtual GCPs for the photogrammetric processing of photos acquired on 8 June 2018 and 13 May 2019 (Fig. 2a). As shown in Figs 2b, c, d, these GCPs are mostly salient features such as bedrock or boulders on the deglaciated banks. Afterwards, the final DSMs and DOMs are within the geographic reference of 17 October 2018. Note that the output resolution is uniformly set to 10 cm per pixel to prescribe a constant resolution of final DSMs and DOMs.

The lack of physical GCPs used to bound the potential distortions of the reference dataset (17 October 2018) would lead to an uncertainty of DSM accuracy. The vertical errors would be further transferred to other products by the virtual GCPs. Therefore, for glacier areas where it is challenging to lay GCPs in situ, UAVs with RTK/PPK can improve the accuracy of DSM and DOM due to their more accurate on-board position and orientation system (POS), hence should be the first choice for glacier surveys.

Uncertainty assessment

Ground checkpoints are usually used to correct the geometry and optimise the accuracy of the 3-D point cloud during the photogrammetric processing, as well as to evaluate the accuracy of photogrammetric aerial triangulation and reconstructed products. However, because of the steep terrain and rough glacier surface in the glaciation zone of the debris-covered Dagongba Glacier, it is challenging to set up and survey dozens of check points during the field campaigns. Ideally physical GCPs, stable areas and other available elevation products should be considered for the uncertainty assessment.

In the stable area, displacements of corresponding points on the multitemporal datasets should be zero in an ideal case, which is one of the most common methods for estimating errors in feature-tracking glacier surface velocity (Kraaijenbrink and others, 2016; Sahu and Gupta, 2019). If the displacement of glacier-free stable areas in the DOMs and DSMs on 8 June 2018 and 13 May 2019 is not zero, then an error should be assumed in the data processing. Considering the effects of seasonal changes on vegetation, selected points in the bare surface area (over 10$^{5}$ pixels) of the glacier-free zone (Fig. 2a) are manually selected as the statistical sample for error estimation to quantify the uncertainty in this study.

Flow velocity

The co-registration of optically sensed images and correlation (COSI-Corr) is a software package integrated into the Environment for Visualizing Images (ENVI), which provides tools to orthorectify and co-register optical remote-sensing images and to perform automated subpixel correlation (Leprince and others, 2007), making it possible to determine ground surface displacements (e.g. coseismic deformation, ice flow, or slow landslide) from multitemporal images (Leprince and others, 2008). Its correlation algorithms can be based on topographic knowledge from images obtained by different platforms (satellite, aerial and ground-based) without the use of external information such as GCPs. Here, we apply COSI-Corr's frequency correlator to high-resolution UAV imagery (DOM) to derive the surface velocities of the glacier.

The glacier tongue is commonly debris-covered and heterogeneous, which leads to significant changes in multitemporal DOMs that are unrelated to the surface displacements in the ablation area, such as debris slope slumping, englacial conduit collapse and ice cliff calving. Therefore, considering temporal decorrelations (drastic changes in a period), shadowing differences, or other kinds of noise, the correlation process should be as insensitive as possible. COSI-Corr's frequency correlator has the potential to improve this problem caused by these unwanted disturbances in its multiscale mode (Ayoub and others, 2009). This procedure obtains horizontal ground displacement by estimating phase differences between the sliding correlation window in the Fourier domain. It follows a twofold process. The first step is to set up a preconfigured initial window of a larger size to expect displacement. Once the pixel-level relative displacement is approximately estimated, the second step is performed to estimate the subpixel-level displacement, which uses a smaller window size (Leprince and others, 2007).

In this research, an initial window size of $512\times 512$ pixels and final window size of $128\times 128$ pixels are the optimal configurations through trials, which would produce good correlation results (Kraaijenbrink and others, 2016). To ensure sufficient output detail, a step size of 16 pixels is employed to work with the 0.1 m grid size data. Note that the output resolution of displacement maps is consequently 1.6 m. Each output correlation result includes two correlation images and a signal-to-noise ratio (SNR) file. The two correlation images represent the relative displacement components of the glacier (north/south and east/west) in DOM. The SNR file is determined by each measurement to assess the correlation quality, ranging from zero (failing correlation) to one (optimal correlation). Afterwards, the displacement maps along the lines estimated from the east/west image and along the columns estimated from the north/south image can be obtained. All displacements are then converted to normalised daily values, and the surface velocity during the warm season is analysed from 8 June to 17 October 2018 (a total of 130 d). The cold season is from 17 October 2018 to 13 May 2019 (a total of 208 d).

Elevation changes

For a certain geographic location, the DSM differences do not compare the same features on the surface at the same geographic positions. The observed results for each grid point corresponding to the glacier surface are therefore a combination of net ablation and displacement due to glacier flow (Watson and others, 2017; Brun and others, 2018):

(1)$$W_{{\rm s}} = U_{{\rm s}}\tan\alpha + W_{{\rm e}}$$

where $\alpha$ is the surface slope, $U_{s}$ refers to the horizontal surface velocity and $W_{s}$ is the vertical velocity. The emergence velocity $W_{e}$ refers to the upward flow of ice relative to the glacier surface (Cuffey and Paterson, 2010; Immerzeel and others, 2014; Brun and others, 2018). To avoid bias caused by ice flux, we introduce the horizontal surface displacement derived by COSI-Corr from the DOMs to evaluate the elevation changes for the flow correction, and estimate the mean emergence velocity and downslope motion of the Dagongba Glacier for June 2018–May 2019 with the flux gate method (Vincent and others, 2016; Brun and others, 2018):

(2)$$\left\{\matrix{x_{t + {\rm d}t} = x_{t} + U_{x}{\rm d}t\cr y_{t + {\rm d}t} = y_{t} + U_{y}{\rm d}t\cr z_{t + {\rm d}t} = z_{t} + W_{s}{\rm d}t }\right.$$

(3)$$U_{{\rm s}}\tan\alpha = {z( x + U_{x}{\rm d}t,\; y + U_{y}{\rm d}t) -z( x,\; y) \over {\rm d}t}$$

where $U_{x}$ and $U_{y}$ are the $x$ (east/west) and $y$ (north/south) components of the surface horizontal velocity, $z$ is the glacier surface elevation and ${\rm d}t$ is the duration between different acquisitions. When ice flows along a longitudinal gradient rather than a rough local surface slope, we extract $z$ from the DSM acquired in October 2018 and use filtering for smoothing to reduce the effect of the rough surface (Brun and others, 2018).

We assume that $W_{{\rm e}}$ is homogeneous over the entire glacier tongue. The mean emergence velocity downstream of a cross section can be calculated as the ratio of the ice flow through the cross section divided by the area of the glacier downstream of the cross section (Brun and others, 2018):

(4)$$W_{{\rm e}} = {\Phi \over A_{{\rm T}}}$$

where $\Phi$ and $A_{T}$ are the ice flux and the glacier area downstream of the cross section. This method requires an estimate of the ice thickness along the cross section to calculate the ice flux across the glacier cross section. Due to the lack of field observations, we used the laminar flow equation to estimate ice thickness (Cuffey and Paterson, 2010; Gantayat and others, 2014):

(5)$$H = \ ^{4} \hskip-7pt \sqrt {{1.5U_{{\rm s}}\over Af^{3}( \rho g \sin\alpha) ^{3}}}$$

where the creep parameter $A = 3.24\times 10^{-24}$ Pa$^{-3}$ s$^{-1}$, the scale factor $f = 1$, the ice density $\rho = 917\;$kg m$^{-3}$ and acceleration of gravity $g = 9.8$ m s$^{-2}$ (Cuffey and Paterson, 2010; Gantayat and others, 2014). The ice flux can be calculated by multiplying the depth-averaged profile velocity ($\bar {U} = 0.8U_{{\rm s}}$) and the area of cross section, with the assumption that the mean surface velocity is 80% of the centreline velocity (Cuffey and Paterson, 2010; Brun and others, 2018).

Ice cliff changes

To calculate the elevation changes of the ice cliffs near the supraglacial lakes from the UAV dataset, we manually interpreted and delineated their outlines by combining the DOMs of June 2018, October 2018 and May 2019. Specifically, we mapped their boundaries for all ice cliffs in each DOM, and then topologically merged the two polygons corresponding to the cold or warm seasons of each ice cliff. Finally, the two-phase ice cliff masks corresponding to the cold season (June 2018 to October 2018) and warm season (October 2018 to May 2019), respectively, were applied to their elevation change calculations.

Historical data

Earlier observations concerning glacier movements in this area were reported by Li and Su (1996) based on a network of ablation stakes drilled in the ice surface as reference points. Six profiles were placed on the ice tongue of the Dagongba Glacier, and observed at the end of May and August (1982–83) to obtain stakes displacement during the warm/cold seasons. These measurement points were obtained by an optical theodolite to perform repeated forward intersection measurements with the ablation stakes buried in the ice surface as the marker points. The glacier surface velocities were derived after calculating the relative displacements of these marker points (part of their results are presented in the supplementary material). We visually interpreted the spatial locations of the marker points at that time. To reduce the bias of point measurement, we manually selected a template window ($7\times 7$ pixels) to calculate the mean rates for comparison with our study results.

It is challenging to determine the robustness of this assessment of the long-term decline in glacier velocity. The ITS_LIVE dataset (Gardner and others, 2019) of glacier surface velocities was used to examine the long-term decline of the Dagongba Glacier. Surface velocities are derived from Landsat 4, 5, 7 and 8 imageries using the auto-RIFT feature tacking processing chain described in Gardner and others (2018). Because the coarse resolution (240 m) is insensitive to seasonal fluctuations of the ITS_LIVE velocity data, it is difficult to compare with the UAV data and the findings of Li and Su (1996) at the same magnitude. Therefore, we extracted annual velocity profiles for the period 1988–2017 along a manually delineated centreline profile of the Dagongba Glacier (from 3950 to 4250 m a.s.l.) to provide an examination of its long-term ice velocity change.

Results

UAV mapping accuracy

Figure 3 shows the horizontal (DOMs) and vertical (DSMs) errors of the UAV products obtained after photogrammetry of the sample statistics of the nonglacial zones marked by Fig. 2. The distribution of the output accuracy shows that, compared with the reference position (October), the mean horizontal displacements are $0.05\pm 0.13$ m (east/west $0.05\pm 0.09$ m and north/south $-0.02\pm 0.10$ m) for the warm season, $0.06\pm 0.19$ m (east/west $0.02\pm 0.13$ m and north/south $-0.06\pm 0.14$ m) for the cold season. For the vertical error, the mean vertical displacements are $-0.65\pm 0.44$ m for the warm season and $0.19\pm 0.45$ m for the cold season. Considering that the two periods are 130 and 208 d long, respectively, the uncertainty in horizontal velocity is $\pm 0.10$ cm d$^{-1}$ in the warm season and $\pm 0.09$ cm d$^{-1}$ in the cold season. For vertical uncertainty, $\pm 0.34$ and $\pm 0.22$ cm d$^{-1}$ in the warm and cold seasons, respectively.

Fig. 3. Box plots of horizontal (north/south and east/west) errors estimated from the sample area (locations indicated in Fig. 2a) on the DOMs during the warm and cold seasons. The boxes indicate the interquartile ranges, the whiskers indicate the quartile to extreme ranges and the thick red lines indicate the medians. The histograms and corresponding fitted normal distribution frequency density curves of vertical errors estimated from the sample data on the DSMs during the warm and cold seasons.

Seasonal surface velocities

Fluctuations of the glacier surface velocity are mainly dominated by changes in temperature, subglacial drainage system, meltwater and pressure, with notable seasonal dynamic variations (Xie and Liu, 2010). Figure 4 shows that notably varied flow velocities between the warm and cold seasons on the glacier tongue. Distinct differences are observed between the warm and cold seasons. As shown in Fig. 4c the surface velocity in the warm season is $\sim$1–2 times higher than that in the cold season. The surface velocity distributions in different views and seasons of the flow field are depicted in Figs 5 and 6. In the transverse view, displacement rates along profiles ‘$T_1$’ and ‘$T_2$’ show the raw data (red and blue points) and the average (red and blue lines) over a $\sim$50-m-wide swath.

Fig. 4. Mean daily surface velocity and flow direction (arrows indicate the direction and their length corresponds to the magnitude of the velocity) obtained by frequency cross-correlation for (a) the warm season (area with bad data is masked) and (b) the cold season, with an initial window of $512\times 512$ pixels and final window of $128\times 128$ pixels. Three black lines indicate locations of transversal profile ‘$T_1$’, ‘$T_2$’ and longitudinal profile ‘L’ in Figs 5 and 6, respectively. Numbered measurement points indicate the locations of Table 2, which are used to compare the surface velocity changes in the Dagongba Glacier between 1982–83 and 2018–19 (Li and Su, 1996). (c) Shows the velocity ratio between the warm and cold seasons.

Fig. 5. Plots showing transverse surface velocities along (a) profiles ‘$T_1$’ and (b) ‘$T_2$’. The raw data (red and blue points) and the average (red and blue lines) over a 50-m-wide swath for both seasons.

Fig. 6. Plot shows the raw data (red and blue points) and the average (red and blue lines) of the longitudinal surface velocities along profile ‘$L$’ for the two seasons, respectively. The decreasing gradient at the transition (the junction between northwestward and westward flow directions) part of the glacier tongue is similar during the two periods (dotted rectangle).

In the longitudinal view, along a central flow line in profile ‘L’, the surface velocity almost linearly decreased down-glacier towards the terminus in both periods, which ranged from 1 to 6.5 cm d$^{-1}$ in the warm season and 0.9 to 4 cm d$^{-1}$ in the cold season at 4070 to 4270 m (Fig. 6). However, a larger longitudinal decline of glacier surface velocity is observed in the warm season than in the cold season. In the cold period, the surface velocities of the Dagongba Glacier vary from $\sim$3–4 cm d$^{-1}$ in the upper (eastern) part of the research area to completely stagnant near the terminus. In the warm season, the peak value in the upper region is $\sim$1.5 times ($\sim$7 cm d$^{-1}$) that in the cold season (Fig. 4). Near the lateral-terminal moraine (4150 m a.s.l.), the broadening basin (ice bed) diverges the flow direction to the boundary and slows it down (Fig. 6) (Xie and Liu, 2010). The decreasing gradient at the transition (the junction between northwestward and westward flow directions) part of the glacier tongue is similar during the two periods: 1.5 cm d$^{-1}$ (cold) and 2 cm d$^{-1}$ (warm).

Table 2 summarises the comparison of some measurement points of surface velocities in the findings of Li and Su (1996) and the corresponding results of our study (locations indicated in Fig. 4). For these nine points, the displacement rates decreased to varying degrees, whether it was during the warm or cold season. In 1982, the maximum surface velocity reached 4.67 cm d$^{-1}$ from 4100 to 4200 m a.s.l. during the warm season. In 2018, the maximum surface velocity was only 1.95 cm d$^{-1}$, which was even lower than the flow velocity during the cold season (2.23 cm d$^{-1}$) 36 years ago. At 4200–4300 m a.s.l., the velocity in the 1980s exceeded 11  cm d$^{-1}$ during the warm season, and currently it is <5 cm d$^{-1}$. According to statistics, the mean decreasing rate of measurement points from 1982 to 2018 is roughly estimated to be 58% in the warm season, and 28% during the cold season. In addition, based on the ITS_LIVE dataset of glacier surface velocities (Fig. 7) the ice flow rates along the centreline of the Dagongba Glacier show continuous slowdown since 1988 with the largest velocity changes of the Dagongba Glacier tongue occurred between 1988 and 2004 (maximum mean decrease ice flow rate 5.5% a$^{-1}$ during 1988 to 1992). However, early products were limited by the absence of data and/or the low radiation quality. Especially, before the 20th century, the annual glacier flow velocities at this site had errors over 1 m a$^{-1}$. To show the uncertainty in these results, we indicated the corresponding error ranges in Fig. 7. Although the velocity changes are not always monotonic, with small accelerations observed in some years and regions, the surface velocity of the Dagongba Glacier tongue has shown an overall decreasing trend in the past few decades. Fewer observations of glacier velocity were used for examination, but comparisons with UAV data and Li and Su (1996) yielded similar results, with a 14% decade$^{-1}$ decrease in our results (1982–2018) and a 20% decade$^{-1}$ decrease in the ITS_LIVE dataset (1988–2017).

Fig. 7. Annual velocity profiles (the resolution is 240 m, and the sampling frequency is per pixel) along the Dagongba Glacier centreline (from 3950 to 4250 m a.s.l.) for the period 1988–2017 (velocity data provided by the NASA MEaSUREs ITS_LIVE project) (Gardner and others, 2019). The error bars are the maximum uncertainty ranges of the ITS_LIVE velocity.

Table 2. Surface velocity changes in the Dagongba Glacier

The 1982–83 velocities were derived from Li and Su (1996). The 2018–19 velocities were derived from the UAV-based maps (locations marked in Fig. 4). $R_{\rm d}$ is the decreasing ice flow rate between 1982–83 and 2018–19.

Seasonal elevation changes

We estimated the mean emergence velocity in the study area using the location ‘$T_2$’ in Fig. 4 as a cross section of the ice flux and obtained a result of $\sim$0.03 cm d$^{-1}$, which can be ignored compared with the rate of melt. In addition, the elevation changes caused by down-slope motion were shown in Fig. S8. The elevation change ranges $\sim$0 to $-1$ m from 4070 to 4230 m a.s.l. Particularly in the upper part of the glacier tongue, the downslope motion influenced by the faster flow brought about an elevation change of more than 1 m.

The spatial distributions of glacier surface elevation with seasonal changes corrected from glacier flow are shown in Fig. 8, and their corresponding histograms are shown in Fig. 9. During the warm season, we observed a significant and widespread decrease in surface elevation in most of the ablation zone. Several areas with notable vertical deformation are mostly in the vicinity of supraglacial lakes and neighbouring ice cliffs (black rectangles in Fig. 8c). During this period, the total net change of all pixels over the glacier area in Fig. 8a is $-15.73\times 10^5\pm 2.47\times 10^5$ m$^3$ with a mean elevation change of $-2.81\pm 0.44$ m. For the cold season, the total net mass accumulation is $-5.78\times 10^5\pm 2.52\times 10^5$ m$^3$ over the glacier area in Fig. 8b, with a mean elevation change of $-0.72\pm 0.45$ m. However, the ice cliffs near the supraglacial lakes still show a distinct elevation decline in the area from 4100 to 4200 m a.s.l. From June 2018 to May 2019, the ice surface elevation change throughout the year is the sum of seasonal changes within the same region. The total net change of all pixels over the glacier area in Fig. 8c is $-19.25\times 10^5\pm 3.53\times 10^5$ m$^3$ with a mean elevation change of $-3.46\pm 0.63$ m. Considering the absence of physical GCPs for the DSM (October 2018) acquisition, the vertical orientation errors would be introduced to other DSMs, which can increase the uncertainty of glacier surface elevation changes. Thus, the undulations in surface elevation changes in Fig. 8 may therefore be partly due to biases caused by elevation distortion and partly related to imbalances in seasonal glacier surface ablation conditions (Yang and others, 2020).

Fig. 8. Ice surface elevation changes from different datasets of DSMs after flow correction. (a) The warm season from June to October 2018. (b) The cold season from October 2018 to May 2019. (c) The period from June 2018 to May 2019 (area with bad data is masked) and the locations of Fig. 10 (black rectangle). The values of colour bars on each panel show the ranges of uncertainty in surface elevation change, with the upper limit at the top and the lower limit at the bottom.

Fig. 9. Histograms of the surface elevation difference after flow correction at the scale of the warm season, cold season and annual period for the research area.

Supraglacial lake and ice cliff dynamics

Landforms in the debris-covered glacier tongue can be classified into three types: ice cliffs, supraglacial lakes and moraines (Sakai and others, 2000). Figure 8 shows that seasonality and spatial variability of elevation changes are observed on the glacier tongue, among which the intensive melting area is mostly in the vicinity of supraglacial lakes. Pond energy absorption has the potential to cause significant surface and internal ablation (Sakai and others, 2000; Benn and others, 2012; Immerzeel and others, 2014; Miles and others, 2016, 2018). Ice cliffs, generally covered by dust or thin debris and located neighbouring these ponds, melt at a much higher rate than those covered by thick moraine (Buri and others, 2016; Brun and others, 2018). To explore this in detail, we selected five specific areas (locations marked in Fig. 8), and Fig. 10 shows the changes in surface features around the selected sites. As shown in Fig. 10, a typically expanding glacial lake is observed on panel A, surrounded by a circular ice cliff. The width expansion process from 75 to 90 m and 100 m is surveyed using DOMs. Meltwater and rainfall are gathered here and promote the expansion of supraglacial ponds. This process is generally accompanied by calving in the surrounding ice cliffs. The accumulated meltwater has a strong melting effect on the ice body under the ponds until the water percolates down through the englacial conduits. However, the low-lying area produced by the above process are likely to be refilled by meltwater and rainfall in the ablation period of the following year. Panels B and D are two typically discharged periodically supraglacial ponds. Similar processes were described by Benn and others (2012) and Immerzeel and others (2014) for the Lirung Glacier.

Fig. 10. Evolution of surface features around selected locations (background images are the DOM acquired in June 2018). From June 2018 to May 2019 (indicated in Fig. 8). The coloured outlines are manually edited cliff masks. The surface elevation change (after flow correction) for each ice cliff is the warm season in the upper part and the cold season in the lower part. The final black dashed picture is a field photograph of ice cliff E taken by UAV on 13 May 2019.

Panel E was once the largest supraglacial lake in the Dagongba Glacier catchment, formed by upstream meltwater dammed by a lateral moraine lobe. It has been reported that its area has grown rapidly since 2008, with an area of more than 70 000 m$^{2}$ during filling in summer (Zhang and others, 2019). Sakai and others (2000) suggested that the roof of the conduit could collapse, resulting in the formation of ice cliffs, which would accelerate the melting of the debris-covered glacier. Our research shows that the elevation of the lake basin decreased by $\sim$10 m from 2018 to 2019, and some ice cliffs exposed after water erosion at the bottom of the lake basin were found using UAV images. In addition, Fig. 4 shows that the horizontal displacement of this area is almost stagnant. These phenomena indicate that there is a mass of dead ice under the thick moraine of the lake basin, and the erosion of meltwater during the ablation season leads to a severe reduction in ice mass. Using the textured 3-D model, panel C is considered the drainage downstream of panel E. The greatest vertical changes are caused by an outburst of the lake in panel E in July 2018 (Zhang and others, 2019). During the field investigation in October 2018, we found that a drainage entrance of $\sim$190 m$^{2}$ at the bottom of panel E was connected to a tunnel under the exposed ice cliff in panel C and a large gap in the boundary of the terminal moraine (Zhang and others, 2021).

Discussion

The possible mechanism for differences in seasonal dynamics and the impact on glacier motion

The seasonal variations in ice velocities indicate the fluctuations in driving stress (Kraaijenbrink and others, 2016; Benn and others, 2017). Generally, slope is a key topographical factor that affects the driving stress. However, the terrain of the Dagongba Glacier tongue is relatively flat with most of the slope <$20^\circ$. The surface velocity is likely not to be significantly affected by the slope of the terrain compared with the basal motion (sliding and/or internal deformation). The comparative importance of basal sliding relative to internal deformation can be determined by the spatial variation of velocities on the transverse profile (Copland and others, 2009). The monthly average temperature in this area from April to October is above 0$^\circ$ with abundant precipitation (Li and Su, 1996). Especially during the period of ablation (June to October), bulk water (monsoon precipitation and meltwater) drainage to the glacier bed through the englacial conduits, playing a lubricating role and providing the facilitated conditions for basal sliding. Our transverse surface velocity profiles over the glacier tongue (Fig. 5) show that there is indeed a difference in the transverse gradient between the warm and cold seasons. For example, the displacement rate along profiles ‘$T_2$’ shows a rapid decrease near the margins on both sides, suggesting that basal sliding is the dominant movement mechanism for this portion in the warm season. In areas where basal sliding is dominant, ice tends to move in a slab-like motion, with transverse velocity patterns that are generally high and relatively constant in the centre flow, while increasing rapidly away from the sides, and this motion is known as plug flow (Copland and others, 2009; Kraaijenbrink and others, 2016). In contrast, the displacement rate along profiles ‘$T_1$’ shows a different pattern, with transverse velocity increasing gradually in a parabolic mode from the margin closed to the centre (Fig. 5a), which is generally associated with more deformation-led flows (Copland and others, 2009). Although basal sliding remains important in the glacier centre, where subglacial water may concentrate and enhance lubrication, the distinct increase in warm season flow occurs only over the upper part of this region.

Our results support the association between glacier surface lowering and the development of ice cliffs and the presence of supraglacial (Immerzeel and others, 2014; Buri and others, 2016; Watson and others, 2017; Neckel and others, 2017). The observed elevation changes in ice cliffs are generally higher in the warm season (mean rates of elevation changes range from 5.4 to 7.1 cm d$^{-1}$) than in the cold season (mean rates of elevation changes range from 0.7 to 2.0 cm d$^{-1}$), revealing an interaction between the dynamics of ice cliffs and the potentially large seasonal expansion and contraction of lakes (Miles and others, 2016; Watson and others, 2017). However, our research still provides a limited understanding of the fine spatial and temporal dynamics of supraglacial lake and ice cliff, which requires further investigation over a longer period of time.

Slowdown of the Dagongba Glacier in recent decades

Climate change is generally one of the most important factors for the variations in glacier mass balance (Xie and Liu, 2010). In particular, mountain glaciers in many reports have been observed to be retreating as a response to global warming (He and others, 2003; Kaser and others, 2006; Zhang and others, 2012; Yao and others, 2012; Yang and others, 2013; Hugonnet and others, 2021). Figure 11 shows the changes in annual mean air temperature and annual precipitation since 1980 recorded by Jiulong Meteorological Station in Sichuan, which is the closest meteorological station ($\sim$50 km) on the west side of Mount Gongga to the Dagongba Glacier. We observed an upward trend in temperature, and the annual mean temperature increased by 1.22$^\circ$C from 1980 to 2019, with a warming rate of $0.31\pm 0.10^\circ$C decade$^{-1}$, which was equivalent to the same period on the Tibetan Plateau (Chen and others, 2015; Liu and others, 2015). During this period, the annual precipitation varied slightly with a linear trend of $\sim$$-7.2$ mm decade$^{-1}$.

Fig. 11. (a) Annual mean air temperature and (b) annual precipitation variation in different years at the Jiulong meteorological station (29$^\circ$00$'$ N, 101$^\circ$30$'$ E, 2983 m a.s.l.).

In this case, rising temperatures are most likely the main driver of glacier thinning and volume loss, and the slight variations in the amount of precipitation could not compensate for the glacier mass loss due to the increase in temperature at Mount Gongga (Braithwaite and Zhang, 2000; Kaser and others, 2006; Zhang and others, 2012; Yang and others, 2013; Zhang and others, 2015; Zhu and others, 2018). Pan and others (2012) found that the mean annual temperature on the western slope increased by $0.24^\circ$C decade$^{-1}$ (1988–2009) according to decades of climate records from three meteorological stations around Mount Gongga, and the glacier retreat rate is faster than most other regions in China. In response to the increase in temperature, the mass change is generally accompanied by ice thinning. The ice mass largely controls the changes in gravitational driving stress, which affects variations in velocity (Heid and Kääb, 2012).

Dehecq and others (2019) found that in 9 of the 11 surveyed regions in High Mountain Asia glaciers showed a sustained slowdown in ice flow concomitant with ice thinning over the period 2000–17. The glaciers in southeastern Tibetan Plateau where the Dagongba Glacier is located, are observed to have one of the greatest slowdown rates (Pan and others, 2012). In monsoonal climates, the increase in air temperature not only decreases the accumulation of the warm season but also leads to a sharp increase in ablation along with lower surface albedo (Fujita and Ageta, 2000). Thus, summer-accumulation-type glaciers are more sensitive to warming than winter types. Other meteorological variables, such as evaporation, humidity, solar radiation and wind speed, have an insignificant effect in comparison with the mean temperature and total precipitation (Fujita and Ageta, 2000). In addition, the potential causes for the discrepancy of slowdown between the warm and cold seasons are twofold. On the one hand, mass loss caused by intensive ablation of glaciers during the summer results in a significant reduction in kinetic energy. On the other hand, warming leads to a decrease in the freezing period for glaciers, and the lubrication period with liquid water on glacier flow is prolonged, causing the decrease in glacier velocity during the cold season to be relatively less obvious than during the warm season (Benn and others, 2012).

In addition, limited by the battery power of the UAV and its maximum flight altitude, our study only focused on the ice tongue area below 4300 m a.s.l. in elevation. We thus did not cover the flow velocity and ablation of the entire glacier in this paper. We will conduct longer time series and larger-scale research in future by combining multisource satellite images, including SAR and optical sensors.

Conclusions

This paper presents an attempt to explore seasonal dynamics in response to climatic fluctuations of the debris-covered Dagongba Glacier by using a low-cost UAV for the first time. The raw pictures from three field campaigns (8 June and 17 October 2018 and 13 May 2019) were used to generate DSMs and DOMs after photogrammetric processing. The glacier surface elevation and velocity changes were estimated by differencing DSMs and by feature-tracking based on DOMs, respectively.

In terms of glacier surface elevation changes, there was a widespread and unevenly distributed decrease in most of the ablation zone from June 2018 to May 2019. Especially in the warm season (June to October), the Dagongba Glacier surface elevation decreased by a mean of $2.81\pm 0.44$ and $0.7\pm 0.45$ m during the cold season (October to May). The intensive melting area was concentrated in the vicinity of supraglacial lakes and ice cliffs. During these periods, several localities with considerable mass losses were generally accompanied by the expansion of glacial ponds. The interaction between the development of these supraglacial lakes and the melting of the Dagongba Glacier forms a critical factor that intensifies the mass losses.

To analyse the Dagongba Glacier surface velocity related to seasonal variation, transversal and longitudinal views of the flow field were depicted in both seasons. The surface velocity in the warm season was $\sim$1–2 times that in the cold season, and the spatial distribution increased almost linearly with surface elevation. The velocity of the glacier tongue (from 4070 to 4270 m a.s.l.) ranged from 1 cm d$^{-1}$ to 6.5 in the warm season, while in the cold season, it varied from 0.9 to 4 cm d$^{-1}$.

Compared with the findings in 1982, glacier surface velocity changes have slowed drastically in recent decades. The maximum decline was observed to be up to $\sim$75% at $\sim$4150 m a.s.l. in the warm season. Moreover, the effects of global warming on the slowing down of the Dagongba Glacier in the warm and cold seasons were different. The mean declining rate from 1982 to 2018 over the warm season was approximately twice as large as over the cold season. Variations in glacier surface velocity are generally recognised as an indicator of the mass balance, and the significant decrease in the velocity changes of the tongue emphasises the negative mass balance of the Dagongba Glacier in the past 36 years.