Validation and application of sequential unmanned aerial vehicle 1 surveys to monitor the kinematics of a rapid rock glacier

. 7 Accurately assessing landform evolution and quantifying rapid environmental changes are gaining importance in the context 8 of monitoring techniques in alpine environments. In the European Alps, glaciers and rock glaciers are among the most 9 characteristic cryospheric components bearing the most prolonged monitoring periods. This study introduces a rigorous 10 procedure to quantify rock glacier kinematics and their associated uncertainty derived from sequential unmanned aerial vehicle 11 (UAV) surveys. High-resolution digital elevation models (DEMs) and orthomosaics are derived from UAV image series 12 combined with structure from motion (SfM) photogrammetry techniques. Multitemporal datasets are employed for measuring 13 spatially continuous rock glacier kinematics using image matching algorithms. This procedure is tested on seven consecutive 14 (from 2016 to 2019) UAV surveys of Tsarmine rock glacier, Valais Alps, Switzerland. The evaluation of superficial 15 displacements was performed with simultaneous in-situ differential global navigation satellite system (GNSS) measurements. 16 During the study period, the rock glacier doubled its overall frontal velocity, from around 5 m yr ˗1 between October 2016 and 17 June 2017 to more than 10 m yr ˗1 between June and September 2019. Using the adequate UAV survey acquisition, processing, 18 and validation steps, we almost achieved the same accuracy as the GNSS-derived velocities. Nevertheless, the proposed 19 monitoring method provides accurate surface velocity fields values, which allow an enhanced description of the current rock 20 glacier dynamics and its surface expression.


SfM workflow 186
Recent improvements in photogrammetric processing capabilities, together with advances in computer vision algorithms, have 187 facilitated the emergence of SfM with multi-view stereo (MVS) workflows. These developments have been capitalised by 188 several open source and commercial SfM software packages, such as MicMac, Agisoft Photoscan, 3DF Zephyr and 189 Pix4DMapper, among others (Smith et al., 2016). This study applied the SfM workflow as implemented in the commercial 190 software Pix4DMapper Pro version 4.4 (https://pix4d.com/pix4dmapper-pro/, last access: 3 February 2020). This software 191 provides a straightforward pipeline processing from raw images acquired by UAV devices to point cloud and orthomosaic 192 generation through mainly three steps (Fig 2a). DEMs using an Inverse Distance Weighted (IDW) interpolator, whereas the ensemble of oriented images is orthorectified and 209 mosaicked to generate a true-distance colour orthomosaics at 0.1 m pixel size. 210

Surface movements from sequential orthomosaics 211
Image matching using a normalised cross-correlation (NCC) function on CIAS software (Heid and Kääb, 2012; Kääb and 212 Vollmer, 2000) was applied to sequential orthomosaics (obtained from the previous step) covering the Tsarmine rock glacier 213 and its environs. This procedure relies on the heterogeneity produced by the shape and size of the boulders existent at the rock 214 glacier surface, providing high-quality contrasting targets. Basically, the NCC matches homologous points from the rock calculated for potential homologous points of the reference window within the search window. The homologous point that 218 obtains the highest correlation value is established to be the new point location, and therefore the 2-D displacement, from time 219 1 to time 2 (Kääb and Vollmer, 2000). 220 221 Horizontal surface displacements covering the Tsarmine rock glacier were derived from consecutive orthomosaics resampled 222 to one-tenth of the original pixel size (i.e. 0.01 m). Resampled images allowed achieving displacements at sub-pixel precisions 223 (Debella-Gilo and Kääb, 2011). Surface points for image matching were regularly spaced within the landform boundaries 224 based on a 10 m sampling grid in an Eulerian framework. We also included the 35 points from each TGS with their 225 corresponding initial coordinates as additional points for sequential CIAS calculations. A reference window (initial 226 orthomosaic) of 128 × 128 pixels and a search window (subsequent orthomosaic) of 256 × 256 pixels were found suitable to 227 compute surface displacements from few centimetres up to around thirteen metres. Before applying CIAS, the seasonal snow 228 cover was masked out due to its interference to accurately track surface boulders. Furthermore, a directional filter was used to 229 all computed vectors, considering the primary orientation and the slope gradient on the rock glacier. Nevertheless, the 230 percentage of spurious mismatches were below 5% of the total matches for all periods. An ordinary kriging interpolation was 231 applied to single data voids (i.e. one or two consecutive mismatches) for filling the data gaps in our data. 232 233 Despite using the same combination of GCPs for the sequential UAV surveys, a coregistration assessment of orthomosaic pairs 234 was performed using a Helmert similarity transform on stable sectors near Tsarmine rock glacier (Fig 1c). The delineation of 235 the stable ground was based on the previous work done by Kummert and Delaloye (2018), who thoroughly assessed three 236 stable sectors with repeated TLS surveys between 2013 and 2016 (their Fig 2). Nevertheless, TGS measurements revealed that 237 the narrow and elongated southern stable sector had moved at about 0.15 m yr -1 between 2016 and 2019. Therefore, this sector 238 was excluded during the coregistration assessment. Estimates of the directional variance and the systematic (bias) errors were 239 calculated for the consecutive orthomosaic pairs using around 69 stable rock surfaces circumscribed to the reassessed stable 240 sectors (Fig 1c). The Helmert similarity transform included the detection of systematic rotations, translation (x-y-shift vector), 241 and scale differences between the orthomosaics. However, we found neither scale nor rotations differences between the 242 consecutive orthomosaics, thus demonstrating the horizontal positional quality of the SfM-derived products. Removing the 243 systematic x-y-shift vector from the CIAS measurements provided a bias-free overall displacement. Similarly to Eq. (2), the 244 uncertainty of the derived horizontal displacements using CIAS can be calculated, taking into account the coregistration's 245 anisotropy in both the x and y directions (Redpath et al., 2013). Furthermore, as the x and y have been corrected for the 246 systematic (bias) coregistration errors, the standard deviation for each displacement measured by CIAS is unique and follow: 247 where and are standard deviations obtained from coregistration analysis of the stable rock surfaces using the Helmert 249 similarity transform. In this study, we multiplied the results from Eqs. (2) and (3) by a factor of 1.645 (i.e. confidence limit of 250 90%), and used them to assess the minimum limit of detection (LoD) for each displacement vector. 251

DEM assessment and elevation change analyses 252
The GCPs and CPs root mean square (RMS) errors from the bundle bloc orientation (Table 2)  whole range of UAV-derived velocities displayed a good agreement with the displacements of these points as measured by 273 GNSS surveys (Fig 3). The correlation values for the consecutive periods were strong and significant, with R 2 values ranging 274 between 0.98 and 0.99, indicating high reliability of the methodological approach. Individual outliers are mostly associated 275 with kinematic points marked near the boulder edge. This is mainly because the cross-correlation tends to mismatch 276 homologous points due to substantial changes in the geometry and lighting conditions. On the contrary, kinematics points 277 installed near the boulder centre provide ideal cross-correlation targets for finding homologous points.  Table 3). The largest errors are associated with the first 281 survey (October 2016). They are mainly due to a lower coregistration quality, which may derive from the strong shadows and 282 a thin sheet of snow existing during the UAV flight (Fig 4a). Still, values below LoD were concentrated at the northern and 283 southern margins (levees), where the rock glacier surface is virtually stable (Fig 4).  (Table 4). On the other hand, the 302 mean surface displacements presented large seasonal differences between the snow cover (October to mid-June) and snow-303 free (mid-June to late-September) periods. During the snow-free periods, displacement values share nearly a quarter of the 304 yearly component of the surface displacement (Table 4). It is worth to mention that these seasonal variations of rock glacier 305 displacements are in line with those observed by in-situ methods (Delaloye and Staub, 2016). 306

High-resolution surface observations 307
Aside from describing detailed surface velocity fields employing high-resolution orthomosaic pairs, the associated UAV-308 derived datasets such as DEMs also allow describing some remarkable geomorphic changes on the rock glacier. The DoD analyses show a strongly heterogeneous elevation change pattern along the rock glacier (Fig 7). Aggregated over the 317 surveyed area, the net elevation change between October 2016 and September 2019 is -0.21 ± 0.07 m. On the lower rock 318 glacier section, large thickness changes varied between -7 to up to 7.5 m. However, almost non-significant changes are 319 encountered in both levees sector. During the same period, the estimated volume gain is 12360 ± 1784 m 3, whereas the volume 320 loss is 22800 ± 2954 m 3 of material. The overall net volume change equates to -10526 ± 3451 m 3 , indicating that wastage of 321 material is considerably larger than the accumulation of material over the same surface. It is important to remark that the upper 322 gully section, where the rock glacier is evacuating sediments from the terminus, is not included in the DoD analyses (for 323 detailed characterization, the reader is referred to Kummert and Delaloye, 2018). 324

325
The 3-year development of the scarp structure across the entire width of the rock glacier (Fig 5), which was already existing 326 at the time of the first survey, can be observed from the dynamic visualisation of the hillshade images at the video supplement 327 (V1). At the time of the last UAV survey, the height of the scarp was reaching locally up to 10m (Fig 8). This dynamic 328 visualisation also shows the passive transfer of material and an almost negligible rotational movement of boulders. 329 Additionally, the video depicts the progressive enlargement in the southwestern margin (Fig 6b), which seats in direct contact 330 with the stable southern levee, whereas in its northwestern portion, the rim has not expanded very significantly. 331 5 Discussions 332

UAV monitoring strategies 333
The close-range remote sensing monitoring approach based on repeated UAV surveys presents some benefits compared to 334 other more classical remote sensing and in-situ techniques. First, as each segment of the remote sensing chain is controlled 335 (Schott, 2007), from image acquisition, processing, and analysis, this approach is highly customisable to different monitoring 336 periods (i.e. temporal resolution) and to the desired level of detail (i.e. spatial resolution). This is not trivial, as the users of 337 classical photogrammetric surveys or satellite imagery cannot directly operate the platform, leaving their management into the 338 hands of commercial or governmental agencies. Furthermore, this may hamper the optimal monitoring period and set 339 prohibitive operational costs for the remote sensing strategy. Second, compared with TGS techniques, the time for data 340 acquisition may be reduced when deploying UAV devices (provided that permanent GCPs are already established). On the a regular TGS may take some hours (depending on the personnel's expertise) to measure the ensemble of kinematic points. 343 With nearly optimal conditions during image acquisition, the LoD from UAV derived velocities is nearly as good as the one 344 from GNSS measured velocities, but it can be substantially degraded during unfavourable conditions. Besides, the use of UAV 345 helps increase personnel security by avoiding the terrestrial surveys on unstable sectors. This is especially the case in Tsarmine, 346 where many boulders have become unstable and where the frontal area is always prone to rockfalls due to the strong rock 347 glacier acceleration. More recently, the work carried out by Fey and Krainer (2020) also employed UAV and GNSS data to derive flow velocities 371 on the Lazaun rock glacier, which is located in the Schnals Valley, Italian Alps. Nonetheless, because some of the work was 372 entrusted to different external companies (i.e. outsourced), neither the UAV nor the GNSS processing steps could be described 373 in their study. Therefore, intrinsic data processing quality for each survey was not available (like our Table 2). From our side, 374 we performed all the data collection, essential processing steps, and analyses of our datasets. Furthermore, they employed the the UAV-derived velocities validation (Fig. 3). In this fashion, we want to emphasise that our kinematic points are explicitly 377 used to provide independent accuracy examinations. Besides, such a large number of GCPs will be impractical for some rock 378 glaciers, where terrain constraints such as unstable ground and accessibility will compromise the safety of the monitoring 379 efforts. 380

UAV and TGS frameworks 381
As we indicated in our methodology using sequential orthomosaics, we consider the calculation of UAV-derived velocities 382 using a fixed grid system (video supplement V2) through which the rock glacier surface flows (i.e. Eulerian specification of 383 the surface flow fields). It should be stressed that during the validation step, we actually employed the initial coordinates of 384 the kinematic points that corresponded to each UAV survey acquisition (Sec  , 2018), and therefore the rock glacier terminus is oscillating from season to season (Fig. 6). During 400 most of the snow cover periods, the freezing of the active layer explains the net frontal advances identified in mid-June. This 401 freezing leads to increasing ice content and cementing rock particles, and consequently, reducing the frontal erosion rates 402 during the cold periods (Kummert et al., 2018). Contrarily, the thawing of the active layer and the freshly exposed permafrost 403 ground during the snow-free periods generates increased erosion rates and net frontal retreats (Fig 6). The high frontal 404 velocities and the enhancement of the transversal scarp feature indicate a destabilisation phase (Marcer et al., 2019). However, 405 a potential collapse of the landform is not expected, even if the high velocities persist, due to two main reasons: (1) the net 406 surface elevation gain observed along the rock glacier front which is leading toward to a concave profile (Fig 7) and; (2) the The onset of the scarp feature predates the monitoring period using UAV surveys, and the use of sparse kinematic points from 410 previous TGS can help reveal the initial conditions of this feature. As can be seen from Figure