Long-range terrestrial laser scanning measurements of annual and intra-annual mass balances for Urumqi Glacier No. 1, eastern Tien Shan, China

The direct glaciological method provides in situ observations of annual or seasonal surface mass balance, but can only be implemented through a succession of intensive in situ measurements of field networks of stakes and snow pits. This has contributed to glacier surface mass-balance measurements being sparse and often discontinuous in the Tien Shan. Nevertheless, long-term glacier mass-balance measurements are the basis for understanding climate–glacier interactions and projecting future water availability for glacierized catchments in the Tien Shan. Riegl VZ®-6000 longrange terrestrial laser scanner (TLS), typically using class 3B laser beams, is exceptionally well suited for repeated glacier mapping, and thus determination of annual and seasonal geodetic mass balance. This paper introduces the applied TLS for monitoring summer and annual surface elevation and geodetic mass changes of Urumqi Glacier No. 1 as well as delineating accurate glacier boundaries for 2 consecutive mass-balance years (2015–2017), and discusses the potential of such technology in glaciological applications. Three-dimensional changes of ice and firn–snow bodies and the corresponding densities were considered for the volumeto-mass conversion. The glacier showed pronounced thinning and mass loss for the four investigated periods; glacierwide geodetic mass balance in the mass-balance year 2015– 2016 was slightly more negative than in 2016–2017. Statistical comparison shows that agreement between the glaciological and geodetic mass balances can be considered satisfactory, indicating that the TLS system yields accurate results and has the potential to monitor remote and inaccessible glacier areas where no glaciological measurements are available as the vertical velocity component of the glacier is negligible. For wide applications of the TLS in glaciology, we should use stable scan positions and in-situ-measured densities of snow–firn to establish volume-to-mass conversion.


Introduction
Glacier meltwater is a crucial freshwater resource for populations and hydro-economies in arid and semi-arid regions (e.g. Sorg et al., 2012;Chen et al., 2016). The concept of 'solid reservoirs' is well represented in the Tien Shan, where most glaciers have experienced substantial mass loss over recent decades (Farinotti et al., 2015;Pieczonka et al., 2015;Liu and liu., 2016;Sakai et al., 2017;Li et al., 2018). Hence, a better understanding of the relationship between Tien Shan glacier wastage and changing climate is 30 important for projecting water availability in the near future. Glacier mass balance provides important information on the gain or loss in glacier mass and is a direct and immediate indicator of climate evolution (Kaser et al., 2006;Haeberli et al., 2007).
Continuous mass-balance observations are fundamental to understand climate-glacier interactions (Zemp et al., 2015). Annual and sometimes seasonal surface mass balance of individual glaciers can be measured using the direct glaciological method. Stakes are 35 drilled into the ice, allowing the monitoring of ablation, and snow pits are dug in the area where snow has accumulated to provide net accumulation (Østrem and Brugman, 1991;Xie and Liu., 1991;Cogley et al., 2011). However, the shortage of long-term financial and human resources and the inaccessibility of remote regions and natural hazards means that ongoing in situ glacier mass-balance measurements are sparse and discontinuous in the Tien Shan, only Tuyuksu glacier (northern Tien Shan, Kazakhatan) and Urumqi Glacier No.1 (eastern Tien Shan, China) have long glaciological mass balance series (Hoelzle et al., 2017). In contrast 5 to the extensive in situ measurement networks required for glaciological observations, the geodetic method provides mass balance by repeated surveys of the glacier surface terrain, in which two digital elevation models (DEMs) are subtracted to calculate the volume changes and then convert them to mass balance using a density conversion (Zemp et al., 2013;Huss, 2013;Andreassen et al., 2016). The method includes all processes that affect the surface, internal and basal mass balances (Cuffey and Paterson, 2010;Sold et al., 2016) but the geodetic mass balances are assumed to be accurate since the topographic surveys are of high quality 10 (Thibert et al., 2008;Huss et al., 2009;Joerg et al., 2012). The available DEMs, derived from aerial photography and traditional remote sensing imagery, usually limit the accuracy and time resolution of geodetic mass-balance measurements (e.g. Cox and March et al., 2004;Cogley, 2009;Fischer, 2011).
In recent years, emerging earth observation technologies (e.g. airborne (ALS) and terrestrial laser scanning (TLS)) allow the 15 derivation of high-resolution DEMs with vertical and horizontal errors on the order of a few centimeters, and increasingly been used to calculate geodetic mass and changes in glacier volume (e.g. Geist et al., 2005;Pellikka and Rees,2009;Joerg et al., 2012;Gabbud et al., 2015;Fischer et al., 2016Klug et al., 2018. ALS is effective for rapidly mapping extensive areas, but the difficulty of studying glacial changing processes with high temporal resolution since the high costs of ALS and the presence of great topographic relief and high-altitude rock outcrops around glaciers reduce the capacity of observations by aircraft as most ALS 20 instruments have limited operating flight altitude, so we need ground-based surveys (Piermattei et al., 2015). The TLS system is usually simpler, more economical and more flexible than ALS, and has become a well-established tool for monitoring annual and sometimes seasonal changes in individual glaciers (e.g. Gabbud et al., 2015;Fischer et al., 2016;López-Moreno et al., 2016). The new high-speed and high-resolution Riegl VZ ® -6000 terrestrial laser scanner offers a long measurement range of more than 6 km and a wide field of 60° vertical and 360° horizontal for topographic (static) applications (RIEGL Laser Measurement Systems, 25 2014a). The scanner is a Laser Class 3B, with laser wavelength in the near-infrared (~1064 nm), and thus well-suited for measuring snow-and ice-covered terrain in repeated glacier mapping. One study has covered the novel use of Riegl VZ ® -6000 TLS to measure surface melt for a temperate Alpine valley glacier at the seasonal and hourly scales (Gabbud et al., 2015); however, only the middle and lower elevations were detected as the glacier is relatively big. Another study reports the performance of the Riegl VZ ® -6000 in monitoring the mass balance of five glaciers in the European Alps; the surface terrain of each glacier can be almost 30 entirely detected using one scan position since these glaciers are very small and have steep terrain (Fischer et al., 2016). For medium-sized and large glaciers with flat terrain, however, one scan position cannot survey the whole glacier surface.
Urumqi Glacier No.1 has the most detailed annual and seasonal surface mass balance measurements in China. It is also one of the reference glaciers in the World Glacier Monitoring Service (WGMS) network due to its long data series, important location and 35 significant local water supply (Li et al., 2011;Zemp et al., 2009). TLS surveys of Urumqi Glacier No.1were initiated on 25 April 2015 for four scan positions (Fig. 1a), and the subsequent surveys correspond to measurement dates for glaciological mass-balance measurements. Multi-temporal high-resolution and high-precision TLS-derived DEMs are therefore available. To date, comparison of glaciological and geodetic mass balances of Urumqi Glacier No.1was reported for the period 1981-2009 at intervals of several years (Wang et al., 2014) and for the period 1981-2015 , but these studies used a series of low-quality topographic 40 maps to calculate sub-decadal and decadal geodetic results. An accurate reanalysis of seasonal and annual glaciological mass balance of Urumqi Glacier No.1 using high-resolution and high-precision DEMs has not been performed. Our previous study has used the TLS to implement two measurements one month apart (25 April-28 May 2015) to get monthly net mass balance of Urumqi Glacier No.1, whereas we simply compared glaciological and TLS-derived geodetic elevation changes of individual stakes, whether agreement between the glaciological and TLS-derived glacier-wide mass balance was pending, potential and shortcomes 5 of such technology applied in seasonal and annual glacier mass-balance measurements in western China had not been discussed; besides we only considered snow/firn densities in the determination of a density conversion, which was used to convert monthly volume change to geodetic mass changes, as an abundance of fresh snow covered the entire glacier surface at the time of the TLS surveys . In fact, the volume-to-mass conversion becomes more challenging over short time periods as meteorological factors change mass balance gradients (Huss, 2013). Several recent studies have used an area-weighting method to 10 calculate the annual density conversion by classifying a glacier surface into bare ice and firn (Fischer et al., 2016;Klug et al., 2018).
But the volume changes in ice and firn/snow usually take place at the same vertical layer for summer-accumulation-type glaciers (accumulation and ablation take part simultaneously in summer months) from our field observations, it is therefore inappropriate for this study to adopt the area-weighting method. Besides, compaction and metamorphosis imply a shift in the vertical firn profile as well as changes in firn thickness and density (Cuffey and Paterson, 2010;Ligtenberg et al., 2011),so assuming no change occurs 15 in the vertical firn density profile over time in the accumulation area is unrealistic (Bader, 1954).
This study takes Urumqi Glacier No.1 as a case and describes the use of the TLS to monitor annual and seasonal geodetic mass balances for two consecutive mass balance years . The aim of this study is thus to established an optimization scheme of volume-to-mass conversion to realize the calculation of TLS-derived geodetic mass changes, to investigate the possible causes of 20 the differences between glaciological and geodetic mass balance. The potential of such long-range TLS to measure mass balance of glaciers in western China is evaluated and several main considerations for a wide application of the TLS in glaciology are put forward.

Study site
Urumqi Glacier No.1is a northeast-orientated small valley glacier, situated on the northern slope of Tianger Summit II (4848 m 25 a.s.l.) in the eastern Tien Shan (43°06′N,86°49′E,Figs. 1a and b). This glacier covered a total area of 1.555 km 2 on 2 September 2015 according to TLS-derived high-resolution DEMs. Intensive glaciological investigations of Urumqi Glacier No.1were implemented in 1959 and then a monitoring station (Tien Shan Glaciological Station) was set up for long-term glaciological massbalance measurements. During the period 1959-2008, this glacier had experienced two accelerated mass loss, commencing in 1985 and 1996 respectively (Li et al., 2011). The glacier separated into two branches in 1993 (Li et al., 2011): the east branch and the 30 west branch of Urumqi Glacier No.1.
Urumqi Glacier No.1is a typical summer-accumulation-type glacier in a continental climate setting (Liu and Han, 1992;Li et al., 2011). The westerly circulation is influenced by the dynamic action of the Tibetan Plateau in the winter months, causing a cold climate with little precipitation in the study site (Han et al., 2006;Huintjes et al., 2010). During the summer month, the Tibetan 35 Plateau becomes a thermal depression and forms a plateau monsoon, which carries warm and humid air from the India Ocean, producing abundant precipitation surrounding the Plateau (Huintjes et al., 2010). These climatic conditions were confirmed by the annual climate records  of Daxigou Meteorological Station, located about 3 km southeast of Urumqi Glacier No.1 at 3539 m a.s.l.; the annual average air temperature was about -5.0°C, and the annual average precipitation was 460 mm. 78% of the annual total precipitation amount occurs from May to August (summer), dominated by solid precipitation (Yue et al., 2017). The climatic conditions decide that the glacier is dominated by weak accumulation from October to March (winter) and the accumulation is higher from April to May; both strong ablation and accumulation mainly take place between June and September (Liu et al., 1997). 5

Principles and key features of Riegl VZ ® -6000 TLS
Riegl VZ ® -6000 TLS is an active laser imaging technique that calculates the distance between the object and the laser transmitter based on time-of-flight measurement with echo digitization and online waveform processing, and consequently the position of the 10 point of interest to be computed (RIEGL Laser Measurement Systems, 2013). The scan mechanism includes a fast-rotating (60-120° from zenith) and more slowly rotating optical head (0-360°). The mirror deflects the laser beam in different directions, thus forming a scan line from consecutive measurements. Meanwhile, the optical head rotates and this scan movement is called a frame scan. A line scan and frame scan generate a view scan using this technique; data collection occurs at a rate of 23 000-222 000 points per second and generates point clouds accordingly (RIEGL Laser Measurement Systems, 2014a). 15 The high-accuracy and high-precision ranging is based on its unique V-line technology of echo digitization and online waveform processing, which allows Riegl VZ ® -6000 TLS to operate even in poor visibility and in demanding multi-target situations caused by dust, haze, rain, snow, etc. (RIEGL Laser Measurement Systems, 2014a).

Terrestrial laser scanning surveys 20
Multi-temporal terrestrial laser scanning data of Urumqi Glacier No.1were collected from four scan positions to achieve maximum coverage, and each scan location was selected from the directions where most glacier surface point clouds would be achieved (i.e. the best possible visibility to glacier surface terrain) (Figs. 1c,d). To avoid ground motion and to obtain accurate coordinates of point clouds, each scan position was fixed using reinforced concrete with a standard GNSS-leveling point. The four scan positions were surveyed using the real-time kinematic (RTK) global positioning system (GNSS, Unistrong E650 instrument) to give the most 25 accurate direct georeferencing and registration. The 3-D coordinates were acquired in the UTM 45N coordinate system in the WGS84 datum. The accuracy of this type of survey has been reported to be within ±1 cm horizontally and ±2 cm vertically, according to previous studies on Urumqi Glacier No.1 (e.g. Wang et al., 2014;Wang et al., 2017).
After the measurements of 3-D coordinates, the Riegl VZ ® -6000 was mounted on a tripod placed in the scan position to survey the 30 glacier surface terrain. The scan parameters and atmospheric conditions are of crucial importance, which directly determined point cloud data quality (point density and coverage) and acquisition time. The laser pulse repetition rate was first set to 50 kHz, then line resolution and frame angle measurement resolution were set to 0.2° to allow a view scan with vertical and horizontal angles were in the range of 60-120° from zenith and 0-360°, respectively. A fine scan is a rectangular field-of-view scan, and the selected field should always cover the entire glacier to guarantee the overlap percentage of four scans was at least 30% (CH/Z 3017-2015(CH/Z 3017- , 35 2015. With each scan, the laser pulse repetition rate was reset to 30 kHz, and the corresponding line and frame resolution were configured as 0.02° to ensure dense points of the glacier surface, except for the scan campaign on 2 September 2015 (Table 1). All scans are performed on sunny days (dry and windless atmosphere) to avoid the influence of precipitation and fog, which can absorb laser pulse and reduce the possible survey distance. Details of the survey parameters are listed in Table 1.

Point cloud processing 5
Raw data were post-processed with RiSCAN PRO ® v 1.81; this includes direct georeferencing, data registration, vacuation and filtering (RIEGL Laser Measurement Systems, 2014b). For all five scan campaigns, four scan positions were used (Table 1).
At first (i.e. direct georeferencing), the TLS data from the different scan positions had to be transformed from the Scanner's Own Coordinate System (SOCS) into a Global Coordinate System (GLCS). The transformation of a point from SOCS into the GLCS was described by Lichti et al (2005) and can be expressed by the following vector equation 10 where ⃗⃗ is the vector of a target in the SOCS; is the vector of the georeferenced target in the GLCS; ⃗⃗⃗⃗ is the vector of SOCS origin in the GLCS, k is the derived azimuth from the scan position to the backsight station and Hence the direct method of georeferencing uses the 3-D coordinates of the scan positions to realize its functions (Lichti et al., 15 2005;Mukupa et al., 2016;Fey and Wichmann, 2017). The accuracy of the method depends on the quality of the measured coordinates. Previous studies stated that the direct georeferencing technique in TLS using global navigation satellite systems (GNSS) is advantageous compare with total stations and the inclination sensors (Paffenholz et al., 2010;Mukupa et al., 2016).
The location of each scan was fixed in the GLCS after direct georeferencing; but the point clouds of the overlapped areas cannot 20 coincide completely due to the influence of orientation. In the second step, multi-station adjustment (MSA) was used for the data registration of each scan position according to the iterative closest point algorithm (Besl and McKay, 1992;Zhang, 1992). When we used MSA, the location of each scan was locked and the orientation of each scan was constantly adjusted in several iterations to compute the best overall fit for them based on least-squares minimization of residuals.

25
Afterwards we combined the overlapped scans in one layer. An octree algorithm was used to the merged layer to produce points with equal spacing to realize point cloud data vacuation (Schnabel and Klein, 2006;Perroy et al., 2010). A terrain filter was then applied to filter out noise and non-ground data due to atmospheric reflections such as dust or moisture, which still occurred despite scanning on fine days (RIEGL Laser Measurement Systems, 2014b). Finally, visual interpretation was also performed to check the data and remove clear visual outliers,, and then glacier surface point clouds with one layer were produced. 30

Geodetic mass balance calculations
As the orientation of each scan was continually adjusted to compute the best fit, the attitude angles of each scan campaign are different. Multi-temporal registration, also called relative registration, set the processed layer of 2 September 2015 as a reference; alignment of other scan campaigns onto the reference layer was finished with iterative closest point algorithms to determine the 35 spatial bias of the multi-temporal scans and extract accurate elevation changes (Revuelto et al., 2014;Gabbud et al., 2015). The relative registered layers were then exported into LAS data format for further processing. Multi-temporal registration of two consecutive campaigns is a crucial step and determines the reliability of TLS-derived surface elevation changes (Revuelto et al., 2014;López-Moreno et al., 2016;Fey and Wichmann, 2017).

5
After the relative registration procedure, interpolation of the processed point cloud data calculated high-resolution DEMs of the study site. The surface elevation change at the individual pixel r was calculated by differencing the TLS-derived multitemporal DEMs with ArcMAP 10.2. The total volume change was determined by summing the elevation change of different time periods, and is expressed as 10 where K is the number of total pixels covering the maximum extent of Urumqi Glacier No.1, and is the pixel size (1 m × 1 m).
The calculated volume change is converted to geodetic mass balance (m water equivalent (w.e. )) following: where ̅ is the mean glacier area of the two acquisition dates t0 and t1, thinking a linear change over time, is the density of 15 water and is the average bulk density (density conversion) of glacier volume change (Thibert et al., 2008;Zemp et al., 2013).

Density conversion
As described above, the geodetic mass balance is calculated based on volume changes, which require a density conversion.
However, the density is difficult to determine; in most studies, it is estimated and not measured. Some researchers assume that no change occurs in the vertical firn density profile over time in the accumulation area and use glacier ice density for the conversion 20 based on Sorge's law (Bader, 1954). Actually, the firn line, firn thickness and firn density all vary, and using the ice density causes an overestimate of mass balance. Huss (2013) recommended a density conversion of 850 ± 60 kg m -3 for the volume-to-mass conversion based on an empirical firn densification model with idealized surface mass balance forcing. But the recommendation is appropriate in the case of a geodetic observation span longer than 5 years; with stable mass balance gradients, volume changes significantly different from zero and a firn area exists. Therefore, several recent studies classify the glacier surface into firn and 25 bare-ice zones and use the volume-weighting method to calculate the annual conversion (e.g. Fischer et al., 2016;Klug et al., 2018).
However, mass balance processes of Urumqi Glacier No.1primarily occur in summer and the glacier is dominated by weak accumulation in winter (Liu et al., 1997). Glacier volume changes in ice and firn usually occur at the same vertical profile according to long-term observations; the surface classification is not applicable in this study. Here we use in situ measured thickness and densities of firn/snow ( ) and ablation stake data (change in ice thickness) to calculate single-point density where = 900 kg m -3 is glacier ice density, and and are the changes in ice and firn/snow thickness, determined from glaciological mass-balance measurements. We extrapolated single-point values to the glacier-wide densities to calculate average bulk density ( ) using the interpolation method (Table 3), and the distributed density conversions of the total glacier were then 35 generated (Figs. 6a, c, e, g).
It is generally true that the density conversion relies on measurements of changes in the firn body, thickness and density of each firn layer being continuous from the top to the bottom of the snow pit, and a stratigraphic description of the firn layers is completed by experienced investigators. Major change processes in the snowpack (e.g. from crystals to grains, descriptive free water content and ice layers, etc.) can be considered in this case (Kaser et al., 2003), and firn compaction assumed to be negligible. We use t volume-weighting method (the weights are the thickness changes of each firn layer and glacier ice) to calculate the firn density 5 ( ) of each snow pit in this study.

In situ measurements
The mass balance of Urumqi Glacier No.1has been observed using stakes or snow pits, since 1959 (Xie and Liu, 2010).
Glaciological measurements broke off during the period 1967-79 and the glaciological data series during this period were 10 reconstructed from correlations with climatic data observed at Daxigou Meteorological Station (Li et al., 2011). The program was re-established using glaciological methods in 1980. No less than 40 ablation stakes were drilled into the glacier and evenly distributed at different elevation bands using a stream drill, despite the fact that the number of stakes has varied from year to year, and snow pits were dug where snow has accumulated (Fig. 1a). The mass-balance year of Urumqi Glacier No.1is defined from previous September 1 to next August 31 (Liu et al., 1997). Usually, from the beginning of May to early September each year, a 15 spatial distribution of single-point ablation (mass loss) or accumulation (mass gain), and snow density (if there is snow cover) were measured by stakes and snow pits at monthly intervals. The net accumulation is measured by digging snow pits at each of the stakes in the area of the glacier where snow has accumulated during the period of investigation; stakes are drilled into the glacier and change in an exposed stake height plus change in snow depth (if snow exists) at two successive dates gives the net ice ablation at this point (Kaser et al., 2003;Xie and Liu, 2010). Hence the measured items include the stake vertical height over the glacier 20 surface, thickness of superimposed ice, and the thickness and density of each snow/firn layer at individual snow pits. Note that fresh snow covered the whole glacier surface at the beginning of the ablation season, so snow pits must also be dug at each of the stakes.

Glaciological surface mass-balance determination 25
Glaciological mass balance includes point and glacier-wide mass balances. The rate of mass gain and loss per unit time is accumulation rate ̇ and ablation rate ̇, respectively, and the different between these two, i.e. ̇ minus ̇ equals mass-balance rate ̇. Integrating ̇ over the time span from t 0 to t 1 gives point mass balance (Cogley et al., 2011): Point values can be extrapolated to glacier-wide mass balance using the contour-line or profile method (Østrem and Brugman, 30 1991;Kaser et al., 2003). Here the time span is often a year or a season, and a seasonal mass balance is classically a winter balance or a summer balance (Cogley et al., 2011). Here t 0 and t 1 are the same as the t defined by Eq. (4). Point values can be extrapolated to glacier-wide specific mass balance using the contour-line or profile method (Østrem and Brugman, 1991;Kaser et al., 2003). For Urumqi Glacier No.1, contour-line and isoline methods had successfully been used to calculate seasonal and annual glacier-wide mass balance (Xie and Liu, 1991), together with simulated values obtained using a simple energy-balance model (the energy divide 35 into shortwave radiation and temperature dependent energy budget) in areas with no measurements (Oerlemans, 2010;WGMS, 2017). This study involves the glaciological measured data over the period 2015-17.

Accurate glacier boundary delineation
Accurate and updated glacier area was important for both geodetic and glaciological mass balance calculations (Zemp et al., 2013).
Fresh snow cover probably led to an overestimate of glacier extent at the beginning of the ablation season; to reduce the influence of snow cover and to extract accurate glacier outlines, we mainly considered glacier extents at the end of the ablation season.
Glacier boundary delineation was performed following Abermann et al. (2010). Firstly three shade reliefs at the end of 5 hydrological years 2015, 2016 and 2017 with an azimuth angle for illumination (300°) were calculated based on multi-temporal high-resolution DEMs to show optimal visualization of contrasts in different aspects. We then delineated the glacier boundary directly by manually digitizing the strongest roughness in the shade reliefs (Fig. 3). Areas of Urumqi Glacier No.1 were 1.555 in 2015, 1.550 km 2 in 2016 and 1.542 km 2 in 2017. Glacier area reduction was primarily attributed to terminus retreat (Fig. 2b, c).

Uncertainties of geodetic mass balance
After multi-temporal registration, errors related to the spatial bias of the multi-temporal DEMs may be negligible. Besides density conversion for converting TLS-derived glacier surface elevation changes to mass balance, uncertainties in the geodetic mass balances derived from the TLS may be related to (1) errors in point cloud data acquisition, including surface terrain and atmospheric conditions (moisture and wind) (Revuelto et al., 2014;Fischer et al., 2016); and (2) errors in data processing and DEM 15 generation, e.g. registration (multi-station adjustment), point cloud vacuation and filtering (smoothing terrain information) (Wheaton et al., 2010;Gabbud et al., 2015;Hartzell et al., 2015).
As mentioned in Sect. 3.1.2, dry and windless days were selected to perform the five scan campaigns. Instability of the TLS influences the registration of single scan positions from each data acquisition campaign, which includes small displacements of 20 scan positions and the vibration of TLS. Each scan position was established on stable rock surfaces using reinforced concrete (the average drilling depth was greater than 80 cm) with a standardized GNSS-leveling point to avoid ground motion. In fieldwork, TLS is mounted using a tribrach on a tripod to level the instrument . Revuelto et al. (2014) found that the vibration of TLS can introduce considerable errors in measurements performed over large scales. In our experience, this issue is mainly relevant to wind, so windless weather conditions are important. Because the registration error cannot be distinguished from the 25 positional uncertainties and the surface, it is difficult to assess registration-induced uncertainty; the error statistics are usually used to evaluate the registration error (Fey and Wichmann, 2017). RiSCAN PRO ® v 1.81 software reports error statistics of the MSA results (RIEGL Laser Measurement Systems, 2014b). The standard deviation of errors ( ) from the set of residuals obtained from registering the point cloud can be considered as an indication of registration quality (Gabbud et al., 2015;Fischer et al., 2016). Values of over stable terrain surrounding Urumqi Glacier No.1are listed in Table 2 for the four periods. Registration 30 quality was higher at seasonal than at annual scales; the higher quality may be attributable to fresh snow cover, which makes the stable terrain smooth.
Despite four scan positions placed at the terminus of Urumqi Glacier No.1, two undetected areas of west branch exist in the (two green polygons in Fig. 5) due to flat terrain ranges from 4050 to 4100 m a.s.l. (Fig. 1a), where the emission laser cannot be 35 received by the laser receiver. We filled these regions using the spatial interpolation method, which can induce potential errors in DEM creations. The lack of dense measured 3-D coordinates of the terrain limits us to assessing terrain-induced errors qualitatively.
For precision, the undetected areas were not taken into account in calculating the geodetic mass balance in fact, related errors were small as the relative proportions of the two areas over the entire glacier surface were minor (3.1% for summer 2015, 3.2% for 2015-16, 3.6% for summer 2016, and 4.6% for 2016-17, Fig. 5). Furthermore, supraglacial river exists at the strong ablation season (June to September due to glacier melting (Fig. 1c, d), which was detected by the TLS surveys. In order to preserve terrain information as much as possible, the octree algorithm built the topological relationship of scattered points to realize the vacuation 5 of the point cloud. Point cloud filtering is also a significant post-processing step because of the dense ablation stake network, which is actually scanned by the device. Fortunately, fine scan generates high-density points of the glacier surface terrain.
No better ways can be used to evaluate the uncertainty of DEMs without precise and well-distributed stable points (Bolch et al., 2017). The standard error ( ̅̅̅̅̅̅̅̅̅) of elevation changes over stable terrain can be considered as a criterion of the uncertainty of the 10 entire glacier (Rolstad et al., 2009;Zemp et al., 2013). The standard deviation ( ) of the stable terrain elevation changes is suitable for estimating the uncertainty of the DEM differences at the individual pixel scale (Fig. 3); in this case the standard error is defined as the standard deviation. However, the spatial auto-correlation must be considered when we calculate the uncertainty of the glacier-wide elevation changes. Thereby, the uncertainty of TLS-derived glacier-wide elevation changes ( ̅̅̅̅̅̅̅̅̅) for individual glaciers were quantified using the geostatistical analysis methods of Rolstad et al. (2009) and written as 15 where denotes the standard deviation of TLS-derived elevation changes over stable terrain. S cor is spatially correlated area.
Given the high density (> 1 point m -2 ) of the TLS data, we can probably assume that the number of independent items is about the number of glacier pixels (cf. Joerg et al., 2012). Here we therefore assume S cor = S. This leads to calculated values of Δ ̅̅̅̅̅̅̅̅̅̅ range from ±0.16 to ±0.25 m (Table 2). 20 Uncertainties related to the density conversion for a single point ( ) were calculated as , where and are uncertainties of ice and firn densities, which were assumed to be ±17 and ±50 kg m -3 , respectively, following Klug et al. (2018). We then extrapolated single-point values to glacier-wide uncertainties ( ) using the interpolation method on the ArcMAP 10.2 platform (Table 3). According to Huss et al. (2009), the uncertainties of the geodetic mass balance 25 ( ) can be estimated using where ̅̅̅̅̅̅̅̅̅ is the average of TLS-derived glacier-wide elevation changes and the related uncertainty relies on the accuracy of the used DEMs.

Uncertainties of glaciological measurements 30
There are additional sources of errors for glaciological measurements, including the harsh climate, which in turn cause uncertainties in glaciological mass balance that are not easy to quantify (Dyurgerov, 2002). Uncertainties in glaciological measurements were classified into three groups: (i) errors in field observations, (ii) errors related to spatial extrapolation over the entire glacier and (iii) errors due to non-updated glacier area. Note that the class (iii) uncertainties appeared to be negligible due to the short time intervals (two consecutive years) in our study. 35 Point measurement uncertainties are prone to errors in stake readings and snow/firn density measurements (Jansson and Pettersson, 2007;Thibert et al., 2008;Huss et al., 2009), sinking or melting-out of stakes and misidentification of the firn layer surface at the end of the last hydrological year (Zemp et al., 2010). Some studies have demonstrated errors of ±0.1 and ±0.3 m w.e. for reading stakes in the ablation and accumulation areas, respectively (Huss et al., 2009). Zemp et al. (2010) determined an overall stochastic uncertainty at ±0.2 m w.e. a -1 for field measurements. Zemp et al. (2013)  The class (ii) errors originate from extrapolating observed values to unmeasured areas, insufficient spatial distribution of measured 15 sites and the interpolation method. Hock and Jensen (1999) evaluated the error of the interpolation method at about ±0.1 m w.e. a -1 for mean specific mass balances. Huss et al. (2009) computed and compared mean specific net balance with randomly reduced annual stake datasets and found that the error was ±0.12 m w.e. a -1 . For Urumqi Glacier No.1, we find that the differences between specific net mass balance at individual sites and in situ measured point mass balance at corresponding sites were in the range of 0-0.04 m w.e. with an average value of 0.01 m w.e., namely, the error of spatial interpolation in the measured area is small. Firn 20 basin and glacier tongue terrain of the west branch are very steep and the upper eastern elevation of the east branch is also precipitous, resulting in no in situ measurements are available in theses inaccessible areas. Therefore the error mainly originates from unmeasured areas (e.g. accumulation areas) but the lack of measured data in the accumulation areas limits us to quantify relevant uncertainties. We conservatively assume that the corresponding uncertainty was ±0.1 m w.e. a -1 (cf. Andreassen et al., 2016). 25 Taking into account the above-mentioned factors, the uncertainty of the glaciological mass balance is calculated as: Resulting values of are listed in Table 3.

Spatial patterns in TLS-derived surface elevation changes
The high-accuracy and -resolution DEMs allowed a detailed insight into the glacier surface elevation changes. Distributed elevation change patterns are generally similar for the four periods, i.e. both branches are characterized by a lower-elevation thinning of 1-3.5 m, elevation changes in the upper-elevation parts are more positive, variations with smaller lowering to pronounced thickening except for west branch in the mass-balance year 2016-17 (Fig. 4a, c, e), this altitudinal change patterns are 35 in good agreement with the long-term glaciological measurements.
Compared to the mass-balance year 2015-16, areas of clearer increase were observed in the upper eastern parts of east branch in the mass-balance year 2016-17, but ice losses in the lower-elevation parts and glacier thickening in the upper reaches of west branch were greater in the first mass balance year (Figs. 5c, g). During summer periods, surface lowering in summer 2015 mainly occurred in the ablation areas of east branch (Fig. 4a), and glacier surface ablation was significantly greater in summer 2016 than in the first summer (Fig. 4e). For a completed mass balance year 2015-16, glacier thinning areas and values in summer were obviously bigger 5 than the whole year, which may be related to fresh snow covered the glacier at the beginning of ablation season. In addations, there were some curves of pronounced glacier surface lowering in the ablation areas during summer periods, which were related to supraglacial river (Fig. 1c, d), a slight thinning area is detected at the lower lift (northerly) edge of east branch, which may be associated with debris cover (Fig. 1c).

Glacier-wide elevation and geodetic mass changes 10
TLS-derived glacier-wide mass balances (Table 3) and their spatial distributions (Fig. 5b, d, f, h) were calculated by multiplying the spatially distributed glacier surface elevation changes (Fig. 4a, c, e, g) with the corresponding distributed density conversion ( Fig. 5a, c, e, g). The thicker snow and firn covered the whole glacier surface at the beginning of the beginning of May each year and the ablation area was bare ice or covered by a thin snow layer at the end of the ablation season according to filed observations (Liu et al., 1997;Xie and Liu, 2010), so the changes in ice and firn/snow thickness are observed during the summer months. 15 However, firn and snow densities are far smaller than glacier ice density, these result in annul single-point density conversion is bigger and the glacier-wide annual density conversions were accordingly higher than the summer ones (Table 3). From the density conversions, we can conclude that the summer geodetic mass balance was highly affected by the snow and firn, and the magnitude of the altitudinal variability in the summer mass balance significantly changed when compared to elevation changes (Fig. 4a, c), whereas the spatial distributed patterns between elevation changes and mass balances showed good performances. These suggest 20 that the density conversions vary for all of the studied periods and a constant value used as the conversion is clearly inappropriate.
Here in situ measured densities from snow pits improved the accuracy of TLS-derived geodetic mass-balance calculations and therefore provided exceptional level of detail on glacier-wide mass balance.
Urumqi Glacier No.1 experienced negative surface elevation changes and mass balances for all of the four investigated periods 25 (Table 3). Summer elevation lowering and mass loss were slightly greater than annual decreases, which may be related to the climatic conditions observed at Daxigou Meteorological Station (see discussion in Sect. 2). In the mass-balance year 2015-16, calculated glacier-wide geodetic mass balance was -0.72 ± 0.17 m w.e., which was slightly more negative than in the second massbalance year. Summer and annual mass balances of west branch were more negative compared to the corresponding values of east branch, except for summer 2016 when the mass loss of east branch was greatest. 30

Comparison to in situ glaciological measurements
TLS-derived geodetic elevation changes at individual stakes closely matched the glaciological elevation change (changes in stake height) of individual stakes from in situ measurements and the difference ( -) in surface elevation changes were close to zero for most of the point measurements (Fig. 4a, c, e, g). The correlation coefficients (R 2 ) between the glaciological elevation change at the ablation stakes and TLS-derived geodetic elevation change at corresponding points were more than 0.90 (Fig. 4b, d, f, 35 h). Note that the location and number of ablation stakes varied slightly over time due to stake melt-out, thus we selected the bestmonitored single-point results to objectively assess the accuracy of the geodetic elevation changes. The majority of the point elevation changes from TLS measurements were slightly positive compared to the glaciological ones, except for summer 2016.
During summer 2016 and mass-balance years 2015-16 and 2016-17, the mean values of -were 0.18, 0.25 and 0.14 m, respectively, which are systematically less than the corresponding uncertainties ( ̅̅̅̅̅̅̅̅̅) of the glacier-wide elevation differences.
The varying tendencies of the glaciological mass balances coincided with the geodetic ones (Table 3). In 2016-17, the difference ( ) in glacier-wide mass balances of Urumqi Glacier No.1between the glaciological and geodetic methods was 5 close to zero. Significant differences between the two methods were detected in summer 2015 for Urumqi Glacier No.1and east branch, with = -0.24 m w.e. and = -0.27 m w.e., respectively. In other three periods, the differences were much less the uncertainties of , which were calculated based on the law of error propagation ( √ ). Overall, the differences of Urumqi Glacier No.1were slightly smaller than those of the two branches (Fig. 6). In order to calculate the statistical significance between the two methods and validate the geodetic against the glaciological mass balance, the reduced discrepancy ( ) between the 10 two methods was calculated following Zemp et al. (2013, Eqs. 19-21); the results of range from -1.14 to 0.58. As falls within the 95% (| | < 1.96) and 90% (| | < 1.64) confidence interval, good agreement between the glaciological and geodetic methods can be considered as satisfactory.
Spatially distributed differences between glacier-wide glaciological and geodetic mass balances were calculated to give the spatial 15 deviations. Over most parts of the glacier surface, especially for the areas near the best-monitored points, deviations were small, indicating both methods showed very close spatial results. Pronounced differences mainly occurred on the steep slopes where in situ measurements were missing (Fig. 7a, d, g, j). The mass balance elevation distribution derived from the two methods remained similar despite the presence of differences in magnitude, i.e. mass balance increased with rising altitude (Fig. 7b, c, e, f, h, i, k, m).
The geodetic results were more positive in lower-elevation regions and more negative in the higher glacier parts in general 20 compared with the glaciological mass balance, which were probably related to glacier dynamic processes (discussed in Sect. 6.4).
The dotted (glaciological) and solid (geodetic) lines met where the glacier mass balances were close to zero; this meant that the equilibrium-line altitudes (ELAs) derived from the two methods matched closely, especially in mass-balance year 2015-16 and summer 2016 (Fig. 7e, f, h, i), but the biggest shift between the two methods was detected in summer 2015 for east branch, which may be related to survey data differences between the glaciological and geodetic observations (see details in Sect. 6.5). This 25 reflects that the TLS can be therefore considered as an effective tool to calculate ELA.

The quality of point cloud data and DEM differencing
The important factors for scanning high-quality point cloud data are visual angles of the scan positions and atmospheric conditions. A dry and windless atmosphere is a prerequisite for high-quality data acquisition. Good visual angles can easily be achieved for 30 very small cirque glaciers. Generally, the area and length of reference glaciers are greater, with a huge variation in altitude. The maximum working distance (6 km) of Riegl VZ ® -6000 is specified for flat targets with size in excess of the laser beam diameter, perpendicular angle of incidence, and atmospheric visibility in excess of 23 km. In bright sunlight the operational range may be considerably shorter than under an overcast sky (RIEGL Laser Measurement Systems, 2014a). However, glaciers generally have complicated surface terrain and the requirement of perpendicular angle of incidence is not always met, so the unscanned regions 35 usually have flat terrain (Fig. 1d). It is very difficult for us to get a dry and windless atmosphere under an overcast sky around a glacier. In these situations, more than one scan position must be set in order to scan as much of the glacier surface area as possible.
However, this, in turn, can create errors in data registration. The average error originating in MSA (σ MSA ) of the investigated periods was ±0.16 m (Table 2). Actually, σ MSA is highly dependent on the overlap percentage of point clouds from each scan position of the same survey data and the accuracy of global 3-D coordinates of each scan position. In our fieldwork, the overlap percentage of point clouds from the four scan positions was more than 30%, which met the requirements of data registration (CH/ Z 5 3017-2015, 2015). Note that we did not get better results with more scan positions since this would probably decrease the quality of the MSA. We should find the best visual angles to obtain the maximum scan range with the fewest scan positions. Higher elevation favors better angles, but it is not always easy for us to access higher parts and place the instrument. For the ideal distribution, scan positions should be located in different elevation bands and directions (Fig. 1a). In addition we can also mount a steadying bar on a tripod to raise the altitude of the TLS. 10 Systematic shifts of DEMs in the horizontal and vertical directions can also increase the uncertainty of DEM differencing (Nuth and Kääb, 2011), so multi-temporal registration of two consecutive scan campaigns is predicated on the TLS-derived geodetic elevation changes being accurate. The mean uncertainty of elevation changes was ±0.22 m (Table 2), which was slightly smaller than in the TLS datasets used for other glacier thinning measurements (e.g. López-Moreno et al., 2016). This may be attributable to 15 the use of accurate global 3-D coordinates of each scan positions and a sufficient number of stable terrains (Fig. 3). In addition, fixed scan positions also reduce multi-temporal registration error and enhance the accuracy of glacier-wide elevation changes. So the quality of point cloud data and DEM differencing is encouraging.

Accuracy of geodetic and glaciological mass balances
It is obvious that the quality of TLS-derived geodetic mass balances relies on the accuracy of glacier surface elevation changes and 20 density conversion of volume to mass changes. With regard to density conversion, our approaches account for the changes in ice and firn/snow thickness as well as the corresponding densities to calculate more accurate values of density (Table 3) (Table 3) and were slightly greater than other geodetic mass balance measured base on the same TLS device (e.g. Fischer et al., 2016), which may be related to multiple scan positions and a larger scanning range, but smaller than those derived from remote sensing imagery (e.g. Holzer et al., 2015;Barandun et al., 2018).

30
Dense spatially measured sites cover the glacier surface (the average density is about 28 stakes km-2 from 2015 to 2017), except for the inaccessible areas (Fig. 1a), to measure the glaciological mass balance. The mean uncertainty of was ±0.12 m w.e. and mainly originated from spatial extrapolation of point measurements. This value is smaller than most recent studies (e.g. Andreassen et al., 2016;Thomson et al., 2017;Klug et al., 2018). This is probably due to relative smaller area and accompanying higher density of point measurements of Urumqi Glacier No.1 than aforementioned studies (Fig. 1a). Thus the TLS device yields accurate 35 geodetic results and the quality of the glaciological mass balances is also very good. Nonetheless, uncertainties in field measurements and interpolation can potentially contribute the deviations between glaciological and geodetic mass balances. The obvious deviations in the spatial distributed differences were found at the firn basin of west branch and the upper right edge elevations of east branch, two unscanned areas also present big deviations (Fig. 7a, d, g, j). Remarkable differences of east branch in summer 2015 induce poorest match between glaciological and geodetic mass balance elevation distribution (Fig. 7c), this can be explained by the lack of well-measured point data (Fig. 7a).

The influence of internal and basal mass balances
The glaciological method cannot measure internal and basal mass balances, but these processes are implicitly captured by the repeated geodetic surveys. We need to provide a rough estimate of internal and basal mass balances of Urumqi Glacier No.1 to 5 detect their contributions to the differences between glaciological and geodetic mass balances.
Urumqi Glacier No.1is a small and cold glacier, basal sliding and bed deformation of the glacier are negligible since ice temperature of the glacier bed has not reached the melting point, so the glacier has low ice velocity and dynamics, and hardly any subglacial water systems (Huang, 1999;Xie and Liu, 2010;. Previous studies have suggested that internal 10 ablation of poly-thermal glaciers is negligible as the ice motion is small (e.g. Albrecht et al., 2000;Zemp et al., 2010). Besides internal melt due to the potential energy of glacier dynamics is negligible since the low reduced dynamics. Thus, internal ablation of Urumqi Glacier No.1 is weak, and mainly comes from the released potential energy of descending water: where denotes annual discharge of flowing water, g is the gravitational acceleration, L f is the latent heat of fusion. Glacier mass 15 loss mainly comes from ablation area and we assume meltwater originates at half the vertical from the glacier terminus to the ELA, so ̅ and are average equilibrium-line altitude (ELA) (4152 m) and the altitude of the glacier terminus (3775 m), respectively, ̅ is the average glacier area between 2015 and 2017. The cumulative measured glacier surface ablation over the two years was used to determine annual discharge and the value of was estimated to be about 1.4×10 9 kg a -1 . A calculation of B pe = -0.005 m w.e. a -1 is made. Internal accumulation results from refreezing percolating water or the freezing of capillary-trapped water 20 in cold snow and firn, which plays a significant role in cold and continental climates ,such as Storglacä ren, Sweden (Zemp et al., 2010) and Castle Creek Glacier, Canada (Beedle et al., 2014), Urumqi Glacier No.1 is a clod, continental glacier in a continental climate setting (Liu and Han, 1992;Li et al., 2011), Urumqi Glacier No.1 is a clod, continental glacier in a continental climate setting (Liu and Han, 1992;Li et al., 2011), measurement and studies on internal accumulation of the is still blank, so we conservatively estimate internal accumulation to be about 4% of the winter mass balance and the resulting value was about 0.01 m 25 w.e. a -1 . Final calculated internal mass balance of Urumqi Glacier No.1was about 0.005 m w.e. a -1 .
Basal ablation is generally attributed to frictional heat of basal sliding and geothermal heating for mountain glaciers (Thibert et al., 2008;Thomson et al., 2017;Galos et al., 2017). Basal ablation caused by frictional heat is very small since hardly any basal sliding of Urumqi Glacier No.1. Here we mainly consider the contribution of basal ablation from geothermal heat (B gt ), which was estimated using 30 , where q = 0.059 W m -2 is the geothermal heat flux (Huang, 1999), t is the mass-balance period; here we primarily consider annual scale and basal ablation was estimated to be about 0.005 m w.e. a -1 . The calculated internal and basal ablation totaled -0.01 m w.e. a -1 .

35
Finally the total value of internal and basal mass balances was close to zero, which is far less than the difference ( ) between the two methods. This suggests that the contribution of annual internal and basal processes is negligible and does not affect the differences between the two methods.

Glacier vertical velocity component
Geodetic measurements of glacier surface elevation changes include glacier surface mass balance and vertical velocity component (Kaser et al., 2003;Geist et al., 2005). Vertical velocity (w s ) is downward (submergence) and makes dynamic thinning of Urumqi Glacier No.1 in the accumulation area, and in the ablation area, vertical velocity is upward (emergence) and makes dynamic thickening of the glacier. This dynamic process results in the general difference between the elevation-distributed mass changes 5 stated above (Fig. 8) Glacier dynamic thinning and thickening can be detected by subtracting the glaciological mass balances from the geodetic ones. 15 For most of the study periods, positive difference values (thickening) dominate in the lower elevations, especially the glacier tongue, and negative difference values (thinning) mainly occur in the higher parts (Fig. 7a, d, g, j). Positive values across the east branch in summer 2015 may be related to different survey dates between the geodetic and glaciological methods. Now applying reciprocal density conversion to the mass balance differences estimate the submergence and emergence velocities. 20 Here we defined the term submergence is negative vertical velocity and emergence is positive vertical velocity. Variation tendency of the estimated velocities at ablation stakes were found to match from in situ measured ones, especially for east branch (Fig. 8).
Relative bigger differences of west branch were detected in the mass balance years 2016-17 (Fig. 8b, d), which may duo to an avalanche in the upper part during the summer 2017. The firn basin terrain of west branch is very steep and is adverse to mass accumulation, which can also be validated in terms of TLS-derived glacier surface elevation changes (Fig. 4g). Thus pronounced 25 misalign of mass balance elevation distribution curves between the two methods occurred. Considering the errors of estimate and in situ measurements, submergence and emergence velocities can be estimated using the TLS-derived DEMs and glaciological mass balance. The difference in mass balance elevation distribution can be largely explained by glacier dynamic thinning at higher elevations and dynamic thickening at lower elevations.
In fact, the vertical velocity of Urumqi Glacier No.1 is small (Fig. 8), we now discuss the errors of glacier surface elevation 30 changes versus dynamic thinning and thickening. Differences in glacier surface elevation changes derived from the TLS and glaciological measurements were close to zero for the vast majority of the ablation stakes, and corresponding errors in the differences were mostly larger than the difference themselves (Fig. 9). Compared with the errors of measurements, dynamic thinning and thickening of the glacier were minor and negligible. So Riegl VZ ® -6000 TLS can be considered as an effective tool to measure the mass balance of Urumqi Glacier No.1 35

Meteorological and glacier surface terrain considerations
Despite glacier thinning and thickening were negligible, other factors such as meteorological conditions and glacier surface terrain may cause mass balance differences between the two methods. Figure 9 shows daily meteorological records provided by Daxigou Meteorological Station from 25 April 2015 to 28 August 2017. Positive temperature and more than 75% of the annual total precipitation amount occurred during the summer months; this probably resulted in summer mass balances that were slightly more 5 negative than annual ones (Table 3). Although the reduced discrepancies between TLS-derived geodetic and glaciological mass balances fall within the 95%, relative bigger difference (∆B) between the two datasets is detectable for east branch, and the difference is caused by survey data differences. The glaciological measurements of east branch were performed five days before the TLS surveys (28 August-2 September 2015). The total precipitation was 67 mm and daily mean temperatures were close to zero (Fig. 9), besides the daily minimum temperatures were all below 0°C during the five days. The climatic conditions may be 10 responsible for the larger mass-balance differences of east branch. In addition, light snow occurred before the last TLS surveys (28 August 2017), possibly resulting in the noticeable increase in the upper eastern parts of east branch (Fig. 5h).
The differences in mass balance between the two methods were possibly related to the effect of glacier surface terrain. The presence of two minimal unscanned areas in TLS surveys is due to the flat terrain of west branch surface (two green polygons in 15 Fig. 5). The geodetic mass-balance calculations did not include these unscanned areas; this cloud potentially increased the difference between the two methods. Furthermore, these undetected regions located in the ablation area and higher wastage than the surroundings were observed according to glaciological measurements. This may imply that the geodetic mass balances of west branch were more positive than the glaciological ones (Table 3) and a discrepancy in mass balance elevation distributions of west branch were observed at 4000-4150 m a.s.l. Nevertheless, the geodetic method is able to cover the majority of the glacier surface 20 and take the terrain characteristics into account, whereas the glaciological measurements cannot capture all the topographic features despite a dense spatial coverage of in situ observations being applied, and what's more, in situ observations are missing in the firn basin and glacier tongue terrain of west branch and eastern elevations of east branch because of the presence of precipitous terrains in these inaccessible regions (green color in Fig. 11a). The eastern elevations of east branch are dominated by the northwest aspect, and the firn basin have aspects from north to northwest (Fig. 11b), aspects are likely to influence the glacier surface albedo and 25 thereby control the surface change patterns (cf. Yue et al., 2017).

Potential of the long-range TLS applied in glacier mass balance monitoring
This study presents the application of multi-temporal Riegl VZ ® -6000 TLS point clouds in mass balance monitoring of Urumqi Glacier No.1. The long-range TLS can provide high-temporal-spatial-resolution and -accuracy DEMs to allow more detailed insight into glacier evolution (e.g. Gabbud et al., 2015). To take advantage of this and provide more-precise glacier surface 30 elevation changes, it is worth remembering that fixed scan positions are highly important between consecutive scans when using our approach. We should also note that not all glaciers in China are as easily accessible as Urumqi Glacier No.1. For many large glaciers, it is not always easy to fix scan positions using reinforced concrete with a standard GNSS-leveling point, but we can mark stable bedrock outcrop as a scan site. Another advantage of this type of TLS is the long scanning range, and such an instrument could allow most of the glacier surface to be scanned from one or several scan positions, especially for remote and inaccessible 35 glacier areas (e.g. crevasses, steep ice, debris cover, etc.). Therefore the instrument provides a quantitative evolution in spatial coverage compared to glaciological in situ measurements, which can be seen as a beneficial complement to glaciological mass balance, particularly for calibrating inaccessible areas. TLS surveys can also provide updated glacier boundary and surface DEMs.
In addition, we can paste several retro-reflective targets (e.g. reflective foils, corner cube reflectors and retro-reflective paintings) to the surface of each stakes and the targets can be easily surveyed and identified since each of them has a high directivity of the reflected laser radiation, then the location of stakes be determined. All of these parameters are favorable for glaciological massbalance calculations. A combination of glaciological and TLS observations may yield optimum results. TLS-derived geodetic results can validate the distributed glacier mass-balance models as the TLS can provide high spatial and temporal resolution 5 measurements, especially in the strong ablation season, the instrument can be used to investigate daily or sub-daily ablation (e.g. Haut Glacier d'Arolla, Switzerland; Gabbud et al., 2015), which can completely meet the requirements of time resolution for glacier mass-balance models.
One drawback of the TLS surveys is the presence of data voids (unscanned areas), even for very small glaciers (e.g. Fischer et al., 10 2016). This is due to limited scanning angle and complex glacial terrain. An emerging low-cost Unmanned Aerial Vehicles (UAV) has the potential to avoid data voids in glaciological monitoring since the good surveying angle of UAV. Immerzeel et al. (2014) showed that UAV combined with a Structure from Motion (SFM) workflow provide a powerful tool for monitoring mass balance and surface velocity of a Himalayan glacier with high spatial accuracy. From our field experiment at Urumqi Glacier No.1, rarefied air and frequent blustery wind around glaciers usually induce the power of UAV were nondurable, and rock outcrops results in 15 difficult operations of such instrument. Hence we mainly consider using UAV to survey unscanned area, integrating of UAV-and TLS-acquired data can provide the whole glacier surface terrain of interest. Other technology such as terrestrial photogrammetry also has the ability to estimate mass balance, and the quality of photogrammetric estimation is similar to the quality of TLS (e.g. Piermattei et al., 2015;Fugazza et al., 2018). However, the reliable of UAV and terrestrial photogrammetry in glacial environments is more dependent on the natural features (i.e. characteristic image objects) of the surveyed surfaces compared with TLS. The cost 20 of TLS is higher than UAV and ground-based photogrammetric surveys.
From our experience, the monitoring tool is potentially applicable to other glaciers provided that these glaciers have small to medium size and relative steep terrain. According to the second Chinese glacier inventory (Guo et al., 2015), ~83% and ~70% of the total number of glaciers in China has an area smaller than1 km 2 and 0.5 km 2 , respectively, and only ~3% with an area lager than 25 5 km 2 . Riegl VZ ® -6000 TLS has been proved successful in monitoring mass balance of very small glaciers (area smaller than 0.5 km 2 ) in the Swiss Alps (Fischer et al., 2016). Hence, the majority number of glaciers can be theoretically measured using the TLS.
Furthermore, if we assume that these glaciers with an area of ≤1.6 km 2 (approximate area of Urumqi Glacier No.1) and a mean surface slope greater than 23.4° (slope of Urumqi Glacier No.1) have a good visibility to be monitored using the TLS, the number proportion of theoretical appropriate glaciers is ~58.5% ,which evenly distributed at different mountains (Fig. 12). However, it is 30 not always easy for us to monitor all of the appropriate glaciers as some of them locate in remote areas (i.e. far away from road), now we have selected select some benchmark glaciers with easily accessible locations for future application of TLS measurements, the TLS system is therefore has huge application potentiality for glacier mass-balance monitoring in China.
Nevertheless, TLS measurements and point cloud data post-processing are challenges for a broader application. One disadvantage 35 of the TLS is that it requires specific knowledge, skills and experience for its use and data processing. Other limitations of the TLS are related to suitable scan positions for obtaining good visual angles of the glacier surface and stable scan positions for multitemporal registration of repeated scans for change detection. In addition, the uncertainties of density conversion still remain at seasonal and annul scales as in situ measured densities of all benchmark sites are difficult to obtain (very sparse glaciers in China have such detailed observations as Urumqi Glacier No.1). The day when relatively smaller amount of snow on the accumulation 40 area and the absence of snow on the ablation area (i.e. snow line is clearly distinguished) should be chosen to perform TLS measurements. We may use a built-in camera of the TLS to create high resolution panorama images of a glacier (RIEGL Laser Measurement Systems, 2014a), then firn/snow and bare ice areas (i.e. snow line) can be determined (e.g. Barandun et al., 2018).
Area-weighting approach can be used to estimate a density because the lack of in situ measured densities makes volume-weighting approach difficult to extensively use. A density assumption over time intervals (≥5 years) based on physical models is also 5 important as most glaciers in northwest China are cold and multi-thermal.

Conclusions
Urumqi Glacier No.1 is one of the reference glaciers in the WGMS network, a representative glacier in Central Asia and the bestmonitored glacier in China. Here, for the first time, we have presented the potential of a novel long-range TLS to monitor summer and annual geodetic mass balances of Urumqi Glacier No.1. The Riegl VZ ® -6000 TLS has long scan range up to 6 km and is 10 exceptionally well suited for measuring snowy and icy terrains in glacier mapping. We use TLS-derived DEMs to calculate summer and annual surface elevation changes and geodetic mass balances of the glacier for two consecutive years (2015-17) as well as to delineate accurate glacier boundaries.
Our analysis suggests that Urumqi Glacier No.1 has experienced pronounced thinning and mass loss for the four investigated 15 periods. Glacier surface elevation lowering and mass loss during the summer were slightly greater than annual values. Glacier-wide geodetic mass balance in the mass-balance year 2015-16 was -0.72 ±0.19 m w.e., which is slightly more negative than in the second mass-balance year. The majority of TLS-derived geodetic elevation changes at individual stakes were slightly positive, but insignificant compared to the glaciological elevation change (changes in exposed stake height) of individual stakes (R 2 ≥ 0.90). The difference in glacier-wide mass balances of Urumqi Glacier No.1 between the two methods was close to zero in 2016-17 but 20 relatively larger differences were detected in summer 2015 for the whole glacier and east branch, which were related to the presence of fresh snow at the time of TLS surveys. Statistical analysis shows that agreement between the glaciological and geodetic methods can be considered as satisfying. Pronounced differences in spatial distributed mass balance mainly occurred at the steep elevations where in situ measurements were missing, which potentially induce the deviations in mass balance elevation distribution.

25
Despite uncertainties inherent in TLS-derived geodetic mass balances, our results show that the TLS device yields reliable results and is therefore well suited to the study of Urumqi Glacier No.1 since the observed vertical velocity component is small. What's more, the TLS can provide accurate and detailed information on glacier area and mass balance changes, its temporal-spatial resolution allows more detailed insight into the glacier's evolution. The greatest strength of the TLS is the long-range scanning which allows most of the glacier surface to be measured, including areas that are inaccessible for in situ measurements. Use of the 30 TLS-based geodetic method will be an important development since it is clearly a beneficial complement to direct glaciological mass balance, particularly for calibrating the unmeasured areas and validating the distributed glacier mass-balance models. A combination of glaciological and TLS observations may yield the optimum results. What's more, the TLS has application potential for glacier mass-balance monitoring in western China as most glaciers (~83% of the total number) have an area smaller than1 km 2 .
For a wide application of the long-range TLS, we can select some benchmark glaciers with easily accessible locations for TLS 35 measurements, but the presence of data voids and snow is still an enormous challenge, the quality of point cloud and DEM differencing and density conversion over short time intervals should be considered.
Data availability. Glaciological mass balance data related to this study are submitted to the WGMS and will be available at website: http://wgms.ch/. TLS point cloud data are available upon request by email to the corresponding author.
Competing interests. The authors declare that no competing interests are present.  . Spatial distribution of TLS-derived glacier surface elevation changes (a, c, e and g); the numbers represent the differences (in m) between the TLS-derived ( ) and glaciological in situ measured ( ) elevation changes at corresponding ablation stakes ( -). Scatter plots of glaciological elevation change against geodetic elevation change at corresponding ablation stakes are presented, and the quality of fittings in terms of R 2 is also presented (b, d, f and h). Black lines are TLS-derived 5 glacier boundary of Urumqi Glacier No.1 and same as the corresponding boundary of Figure 3. White areas indicate outliers, which we have deleted. Two green polygons indicate areas that have not been detected by the TLS.  . Spatial distributed difference derived from TLS-derived geodetic mass balance minus glaciological mass balance (a, d, g and j), black dots represent the location of well-measured ablation stakes, which are same as Figure 5. The hypsometry (50 m altitudinal ranges) and the glaciological (dotted line) and geodetic (solid line) mass balance elevation distribution for the whole study period; both summer and annual mass balances are shown. Gray horizontal bars indicate the area-elevation distribution of 5 Urumqi Glacier No.1. Note that the spatial resolution of glacier-wide geodetic mass balance in summer 2015 and 2016 was downscaled to 5 m, in mass-balance years 2015-16 and 2016-17 was down-scaled to 7 m to coincide with glaciological glacier-wide mass balance. the letters represent ablation stakes. Note than the summer periods and stakes in the higher elevations were not selected for comparisons due to snow cover reduced the quality of in situ measured vertical velocity.
5 Figure 9. Changing differences between TLS-derived ( ) and glaciological ( ) annual and summer surface elevation 10 changes at individual stakes versus normalized glacier elevation range of east branch and west for the four investigated periods. Note than dash vertical lines represents the uncertainty ranges (√ ) of the changing differences ( ), and grey quadrants indicate theoretical areas for Urumqi Glacier No.1 in equilibrium.     Table 3. Glacier-wide mean of density conversion (ρ) and its uncertainty (σ ρ ) (in kg m -3 ) as well as TLS-derived glacier surface elevation changes ( ̅̅̅̅̅̅̅̅̅ ) (in m). TLS-derived geodetic ( ) and in situ measured glaciological ( ) net mass balance at winter and annual scales are listed (in m w.e.).