Multiple Close-Range Geomatic Techniques for the Kinematic Study of the La Paúl Rock Glacier, Southern Pyrenees

Abstract

Rock glaciers are one of the most representative elements of mountain permafrost. Their study can contribute to modelling climate change and its effect on natural and anthropogenic environments. Therefore, it is crucial to understand the evolution and quantify the changes in these periglacial landforms at a global level. This study aims to present the monitoring work carried out on the Pyrenean rock glacier of La Paúl (42°39′40″N, 0°26′34″E) from 2013 to 2020, employing in situ geomatics techniques to determine the landform surface kinematics accurately. For this purpose, global navigation satellite systems (GNSS), terrestrial laser scanners (TLS), and unmanned aerial vehicles (UAV) photogrammetry techniques were used simultaneously to evaluate their compatibility in quantifying displacements. Based on 2D and 3D analyses, the results demonstrate the high surface activity of the rock glacier, with mean variations reaching 36 cm/year (GNSS) and a distribution of deformations that, although intensified on its western side, are present on the entire surface of La Paúl. This study uses state-of-the-art geomatics techniques to present dependable and updated quantitative data on a periglacial landform’s recent development in under-researched areas, such as the Pyrenean temperate high mountain.

1. Introduction

The study of rock glaciers involves analysing one of the most evident landforms of mountain permafrost [1]. They are periglacial landforms of particular interest to the scientific community and are present in many mountainous regions worldwide (e.g., [2,3,4,5]). The distribution and dynamics of rock glaciers provide information on the permafrost behaviour and the climatic and paleoclimatic evolution of mountain regions [6,7,8]. Hence, they are valuable indicators contributing to understanding the high-altitude cryosphere under climate-change conditions. Such is their relevance that the velocity of rock glaciers is being promoted as an Essential Climate Variable by the International Permafrost Association (IPA)—Action Group ‘Rock glacier inventories and kinematics’, as are annual temperatures and precipitation [9,10]. It is, therefore, crucial to have a detailed and accurate knowledge of the kinematics of these landforms.

The Pyrenees harbour valuable information for studying the cryosphere in temperate environments, especially permafrost [11,12]. Studies on the displacement of rock glaciers in the Pyrenees began in the 1990s. Nevertheless, the number of rock glaciers for which precise quantitative data are available remains limited in this mountain range [11,12,13,14,15,16,17,18].

In this context, the present study examines the kinematics (i.e., three-dimensional changes) of the Pyrenean rock glacier of La Paúl over 7 years of annual monitoring from 2013 to 2020. The spatial data were obtained using various in situ geomatics techniques, with detailed processing and analysis workflows employed to determine the creeping processes precisely. This study presents a preliminary analysis of the compatibility of state-of-the-art in situ geomatics techniques in monitoring rock glaciers [19,20] and current quantitative information on permafrost-related landforms in one of the southernmost regions of Europe.

2. Study Site

The Posets massif is the second highest in the Pyrenees and culminates at the peak of the same name at 3375 m a.s.l. On the northwest face of the Posets peak is located the rock glacier of La Paúl (Figure 1; 42°39′40″N, 0°26′34″E). Composed of frost-shattered debris, its tongue, about 400 m long, extends between 2800 m and 3000 m a.s.l. [12].

The La Paúl rock glacier exhibits resilience in both ideal and marginal conditions in the temperate high mountains. It was dated as pre-Little Ice Age (LIA), as advances of the homonymous glacier (currently about 350 m west of the rock glacier) deformed its northwest side during this period. Flow deformation features, arches and furrows predate the LIA [12]. However, rock glacier dynamics persisted after glacial erosion during the LIA and the subsequent period of retreat, moving downslope to the present day.

Although the La Paúl has been less studied than other rock glaciers nearby (e.g., the Posets rock glacier in [14,21,22]), different studies have shown a frozen layer from the ground thermal regime, bottom temperature of snow cover measurements and geoelectric soundings [12,21,23]. The application of GNSS-RTK techniques has supported this deduction by showing annual decimetric displacements on its surface [18,24].

3. Materials and Methods

The kinematic analysis of the La Paúl rock glacier was carried out using Global Navigation Satellite System (GNSS) receivers, terrestrial laser scanner (TLS), and aerial photogrammetry by unmanned aerial vehicles (AP-UAV) for the acquisition of 3D information of the rock surface over eight consecutive years (

Table 1. Geomatics techniques applied to the monitoring of La Paúl rock glacier.

	2013	2014	2015	2016	2017	2018	2019	2020
GNSS
TLS
AP-UAV

).

3.1. GNSS Surveys

The GNSS surveys were conducted with two Leica GPS1200 receivers and RTK (Real Time Kinematic) positioning methods. A base receiver stationed in each campaign on a permanent reference survey pin would send real-time corrections to the GNSS rover receiver to obtain precise relative coordinates. Simultaneously, the base recorded satellite observations during time intervals up to 8:30 h (2020), which would be post-processed with the aid of the network of permanent GNSS stations of the Spanish National Geodetic Reference Network [25] and the Leica Geo Office v8.4 software to obtain accurate coordinates (ETRS89 reference system, UTM zone 30N projection). The annual use of the reference pin and the post-processing of the data made it possible to define the stability of the topographic base and verify the coordinate system.

The surface control of the rock glacier by GNSS receivers was delimited to the control of 20 metal rods and natural points distributed over the entire surface of the tongue in 2013 (Figure 2). The coordinates of these singular points were subtracted to produce planimetric displacement vectors and altimetric deformation models. The latter consisted of rasters generated from TIN interpolation (Triangulated Irregular Network; [26]) of Z coordinate differences using QGIS v3.18 [27].

3.2. TLS Surveys

The scanner used for the massive data capture on the landform between 2016 and 2020 was the Faro Focus3D X330. It is a portable pulsed laser device with Shift-Phase technology with a range reaching 330 m. The scans were configured at a spatial resolution of 6 mm at 10 m for a measurement speed of about 500k pts/s.

The TLS surveys began with a selection of scan locations around the tongue to capture the maximum surface extent of the rock glacier from 360° sweeps. The locations were similar between campaigns, ensuring a high overlap between scans point clouds (minimum common area measured between two scans equal to 20%). Typically, the same number of scans were distributed along the eastern and western sides of the debris tongue, with an intermediate scan along the longitudinal axis of the rock glacier to improve the overlap between the slopes (except in 2016 and 2017 due to time constraints). This led to the acquisition of up to 10 scans per year (7 in 2016, 8 in 2017, 9 in 2018, and 10 in 2019 and 2020) in the La Paúl fieldwork.

The reference frame for georeferencing the surveys consisted of 1 × 1 m red and yellow cloth targets that were easy to transport and identify. In all the scans performed, at least three targets were placed around the scanner so that one of the targets was visible from two consecutive scans on the same side of the tongue (Figure 3). The aim was to strengthen the geometry of the survey. Nonetheless, some control points were not used in the georeferencing process due to their deficient distribution over the tongue, the reduced visibility of their centre from the scanner or the high error in their coordinates. That resulted in campaigns where between 40% (2019) and 80% (2016) of the placed targets were used in the georeferencing process. The centre of these targets was coordinated using the same methodology employed for geolocating the rods with GNSS-RTK methods.

The resulting point clouds after “Cloud to Cloud” registration and “Target based” georeferencing of the scans with FARO SCENE software (v2019 and v2020) were compared using 3D methodologies to avoid the influence of heterogeneous sampling on the quantification of deformations [28]. The algorithm used for the comparison was the M3C2 (multiscale model-to-model cloud comparison; [29]). This algorithm is specially designed for the computation of distances between TLS point clouds in complex geometry scenarios, and it allows the determination of differences with an associated spatially variable confidence level of 95% (LoD95%). The parameters used in the calculation of the deformations were 1 m diameter for the projection scale of a core point cloud of maximum resolution equal to 5 cm. Normals were determined in multiscale mode, ranging from 4 m to 16 m with a step of 4 m. The registration error was derived from the quadratic sum of the registration and georeferencing errors of the scans.

3.3. AP-UAV Surveys

The photogrammetric survey was performed with the DJI Mavic 2 Pro multirotor UAV (Figure 4), equipped with a 20 MP Hasselblad L1D-20c camera. Both the planning and execution phases of the flights were performed according to the methodology presented in [28]. That involved the use of UgCS PRO v.3.4 software to obtain raw nadir and oblique images (25° to vertical) with an approximate ground sample distance of 2 cm/pix and 2.3 cm/pix, respectively, in a flight path composed of strips parallel and perpendicular to the slope direction of the surface of interest.

The UAV imagery capture workflow was similar in the 2019 and 2020 campaigns, where variations in the number of images (428 in 2019 and 407 in 2020) and area photographed were linked to the management of the flight plan in the field (e.g., small readjustments in the flight path or resumption of flights after battery changes).

The scaling and georeferencing of the photogrammetric models were based on the distribution of 12 cloth targets, like those used in the TLS surveys, on the rock glacier surface. They were visible from the UAV and were located on the edges and in the central parts of the area of interest as long as the physical limitations of the environment allowed it (Figure 4). GNSS-RTK devices coordinated the centre of the targets according to the methodology described previously.

The photogrammetric processing was also based on the workflow described by [28], with the support of Agisoft Metashape Pro v1.5 software. In summary, the pipeline followed (I) the search for homologous points between images (alignment) for the determination of internal orientation parameters (focal length, location of the main point of the photograph, and lens distortion coefficients) and external orientation parameters (location and relative position between cameras); (II) definition of targets with known coordinates in the photographs; (III) removal of noise of the tie points generated after alignment; (IV) densification of the tie-point cloud to obtain the final point cloud; and (V) generation of digital elevation models (DEM).

The resulting 3D point clouds and DEMs were compared using different methodologies. The M3C2 algorithm was used to compare the point clouds, where the registration error was obtained from the quadratic sum of the 3D RMSE (Root Mean Square Error) of the targets in each epoch. The DEMs were compared by subtracting their values and identifying their significant vertical variations using Equation (1) in the QGIS v3.18 environment.

$${LOD}_{95\%}=t\sqrt{({{\sigma }_{Z1}}^2+{{\sigma }_{Z1}}^2)}$$

(1)

where t = 1.96, and ${\sigma }_Z$ the vertical standard deviation of the targets.

4. Results

4.1. Surveys Quality

The quality of the GNSS-RTK measurements to coordinate both singular points on the rock glacier surface (i.e., rods and natural points) and the targets supporting the TLS and AP-UAV surveys was set at ±1 cm in XY and ±1.5 cm in Z on average after post-processing. Therefore, a prudent 3D accuracy of ±2 cm was defined for the measurements made on the rock glacier with this equipment.

In the case of the TLS point clouds, the quality was associated with 3D registration errors, equal to or less than ±24 mm (±20 mm mean) in the five campaigns, and georeferencing errors, with values equal to or less than ±32 mm (±27 mm mean).

The overall quality of the photogrammetric surveys was evaluated using the errors provided by Agisoft Metashape, specifically from the RMSE of the control points used to scale and georeference the surveys. Linked to the quality of the GNSS-RTK measurements on the targets, the 3D RMSE was ±40 cm in 2019 and ±56 cm in 2020, with errors between ±32 cm (2019) and ±55 cm (2020) in the Z component.

4.2. Rock Glacier Surveys

The measurement of well-identified points over the entire surface of the rock glacier with GNSS-RTK devices made it possible to obtain information on a large part of its surface (1.5 ha), although at a reduced resolution (1 point per 1000 m²). The measurement of 20 rods and natural points was regular during the eight field campaigns, with exceptions derived from significant errors in the coordinate measurements.

The TLS and AP-UAV techniques significantly increased the number of points measured on the surface over the GNSS-RTK measurements (Figure 5). Although the scanned surface extent was similar to that covered by GNSS-measured points, the point density in the TLS surveys exceeded 1k pts/m² in the tongue. That was higher than the point density of the AP-UAV surveys (about 655 pts/m² on average). Nevertheless, the latter had a larger surveyed area by capturing the entire surface of the rock glacier and part of its surroundings (19 ha).

The complex geometry of the rock glacier surface significantly affected the surveys. In the case of TLS, occlusions between the laser beam and the surface to be measured presented gaps of different dimensions despite multiple scans (Figure 5). The aerial point of view and the capture of photographs in regular intervals through the UAV flight plans, facilitated a homogeneous sampling of the surface with photogrammetric techniques despite the abrupt geometry of the landform (Figure 5). This methodology even made identifying the glacier rocks in the point cloud possible.

The DEMs derived from the aerial photographs showed similar characteristics to the sampling density and extent of the point clouds. They did not present spaces without information or artefacts after interpolation and generation of the elevation models (Figure 6). Like in the point cloud, it was possible to identify blocks throughout the rock glacier surface in DEMs of resolutions of 5.5 cm/pix (2019) and 4.2 cm/pix (2020).

4.3. Rock Glacier Changes

The surface displacements of the points monitored by GNSS-RTK presented mean values of 33 cm/y in planimetry and −20 cm/y in altimetry (

Table 2. Planimetric and altimetric displacements of the GNSS monitored points in the La Paúl rock glacier.

	2013–2020
Point ID	Velocity (cm/y)		Total Displacement (cm)
	XY	Z	XY	Z
1	32	−15	222	−103
2	40	−18	280	−126
3	66	−21	459	−144
4	28	−16	193	−111
5	45	−28	312	−198
6	32	−14	223	−96
7	16	−10	112	−69
8	35	−19	173	−97
9	35	−14	173	−85
10	37	−18	187	−105
11	37	−18	256	−125
12	17	−6	122	−43
13	13	−11	89	−75
14	30	−18	209	−124
15	37	−27	260	−188
16	38	−27	265	−186
17	33	−29	234	−203
18	38	−31	264	−220
19	22	−23	156	−160
20	31	−30	214	−209
Mean	33	−20	220	−133

), causing a total mean deformation in the rock glacier surface of 2.2 m in XY and −1.33 m in Z between 2013 and 2020 (Table 2). There were years with particularly high activity, such as the horizontal downslope movements of 2015/16 (40 cm/y), 2014/15 (35 cm/y) or 2016/17, and 2019/20 (33 cm/y), as well as in the vertical displacements of 2016/17 (−30 cm/y), 2019/20 (−29 cm/y) or 2013/2014, and 2015/2016 (21 cm/y) (Supplementary Materials).

Horizontal displacements were significant in La Paúl (Table 2), reaching maximum mean velocities of 66 cm/y (point 3) and minimum mean velocities of 13 cm/y (point 13). Meanwhile, the vertical component presented more moderate values (Table 2), with maximums of −31 cm/y (point 18) and minimums of −6 cm/y (point 12). The largest 3D displacements occurred at the root (average velocities around 45 cm/y) with decreased velocity toward the front, with displacements around 31 cm/y in the central portion and 35 cm/y at the front (Table 2; Figure 7).

The surface deformations identified from the higher-resolution surveys showed significant displacements along the entire rock tongue. The 3D comparison of the TLS surveys showed annual mean variations of 20 cm/y downslope (

Table 3. Annual 3D variations determined from TLS surveys on the La Paúl rock glacier.

TLS Point Cloud	2016–2017	2017–2018	2018–2019	2019–2020
M3C2 distance (cm; mean ± SD)	20 ± 11	18 ± 16	22 ± 18	18 ± 10

; Figure 8). Deformation magnitude was similar in the different sectors of the rock glacier, with variations of 20 ± 11 cm/y at the front, 20 ± 13 cm/y in the centre and 19 ± 15 cm/y in the uppermost areas. These significant M3C2 distances (LOD95%) were located throughout the surface studied, although they were especially present on the longitudinal axis of the tongue and its western side (Figure 8).

The DEM comparison of the AP-UAV surveys between 2019 and 2020 revealed altimetric losses of −17 cm/y (LOD95%) and deformation patterns similar to those obtained in the M3C2 comparison of the AP-UAV point clouds (Figure 8). On average, significant 3D deformations had a magnitude of 22 ± 6 cm/y (2019–2020), with values ranging from 25 cm/y at the front, 21 cm/y at the centre and 19 cm/y at the root of the rock glacier.

5. Discussion

The application of multiple geomatic techniques at La Paúl between 2013 and 2020 has made it possible to determine the significant surface activity of the rock glacier. The findings indicate that the rock glacier exhibits a creep behaviour characterized by significant horizontal displacements (Table 2). The flow directions of the singular points denote that the creep is not oriented toward the actual rock glacier front (Figure 7) but toward the eroded side in the northwestern portions (where the disequilibrium is greater due to LIA erosion). The tongue thinning is more prominent on the longitudinal axis and western side of the rock glacier (Figure 7 and Figure 8), the flank eroded by the La Paúl glacier during the LIA. These movements were reduced on the eastern side of the rock glacier, with lower velocities.

GNSS-RTK techniques showed velocities of 39 cm/y (2013/20; Table 2), TLS techniques 20 cm/y (2016/20; Table 3) and AP-UAV techniques 22 cm/y (2019/20). According to the IPA—Action Group ‘Rock glacier inventories and kinematics’ [30], these values allow classifying the periglacial landform as active, with average downslope displacements greater than 10 cm/y over its entire surface (Table 2 and Table 3; Figure 7 and Figure 8).

The deformation distribution was consistent between techniques (Figure 7 and Figure 8), as has been shown in other kinematic studies of rock glaciers [19,20]. Nevertheless, the deformation magnitudes showed differences. The −17 cm/y in the Z coordinate obtained from AP-UAV DEMs subtraction between 2019 and 2020 was significantly lower than the −29 cm/y vertical GNSS-RTK displacement (Supplementary Materials). These vertical velocities caused a disparity between GNSS-RTK and AP-UAV monitoring methodologies, although previous studies have shown similar values between techniques [20]. TLS-derived DEMs would not support the Z-coordinate variation due to multiple missing data areas and the abrupt geometry of the rock surface (Figure 6) [28]. Instead, the M3C2 comparison between TLS surveys identified 3D displacements of 20 cm/y in 4 years (2016–2020; Table 3). Values below the 38 cm/y 3D displacement of the GNSS-RTK points in the same period (Supplementary Materials). In addition, the 3D deformations (LOD95%) observed between 2019 and 2020 in the TLS point clouds (18 cm/y; Table 3) and the AP-UAV point clouds (22 cm/y) turned out to be closer, suggesting a better agreement between techniques despite the ground and airborne point of view.

The differences between the deformation magnitudes were associated with the number and location of the measurements obtained by each technique, as previously reported [31]. The 20 points recorded by GNSS receivers in the vicinity of the longitudinal axis differed from the millions of points measured by TLS and AP-UAV techniques over the entire tongue (Figure 5). Despite applying a LOD95%, the high number of measurements, including fewer active zones, may have moderated the average values of surface deformation at La Paúl. In any case, the kinematic analysis showed 2D and 3D displacements in specific zones of the rock glacier concordant between techniques, with values higher than 40 cm/y (Figure 7 and Figure 8).

6. Conclusions

The work carried out on the La Paúl rock glacier over eight campaigns shows an evolution in measurement strategies (i.e., GNSS-RTK, TLS and UAV photogrammetry) and data analysis (i.e., 2D and 3D comparisons). Traditional in situ monitoring techniques and methodologies have been applied and combined with more recent technologies and algorithms, allowing the detailed evolution of the rock glacier surface to be determined at temporal and spatial resolutions otherwise unattainable.

The active rock glacier has shown significant annual displacements, where the horizontal variations on a large part of its surface indicate the creeping processes of its surface at present.

Photogrammetric flights have shown higher cost-effectiveness in working time, surface extension, resolution, and accuracy than other techniques. Future AP-UAV surveys will allow further detailed analysis of the kinematics of the entire surface of La Paúl, which, together with a geomorphological analysis and the application of procedures such as automatic image matching, will improve the knowledge of this landform linked to mountain permafrost in temperate mountains in southern Europe.