Active and inactive Andean rock glacier geophysical signatures by comparing 2D joint inversion routines of electrical resistivity and refraction seismic tomography

In semi-arid Chile, rock glaciers cover a surface area that is four-times larger than that occupied by glaciers. For this reason, their role in freshwater production, transfer and storage is likely to be of primary importance, especially in this area of increasing human pressure and high rainfall variability. To understand their hydrological role now and in the future it is necessary to characterize their internal structure (e.g., internal boundaries, ice, air, water and rock content). In this paper, we present the results and interpretations of electrical resistivity and refraction seismic tomography profiles on an active (El 5 Ternero) and inactive (El Jote) rock glacier in the Chilean Andes. These are the first in situ measurements in Estero Derecho: a natural reserve at the headwaters of the Elqui River, where the two rock glaciers are located. Within our study, we highlight the strong differences in the geophysical responses between active and inactive rock glaciers through the analysis and comparison of three different inversion schemes: individual dataset inversion, structural and petrophysical joint inversion. Moreover, we propose a diagnostic model representation for the differentiation between active and inactive rock glaciers. 10


Introduction
High mountain environments are delicate geosystems under increasing human pressure and represent important geoecological indicators of a changing climate of the areas where they are situated. In particular, in semi-arid Chile (between 29°and 34°S), rock glaciers cover a surface area that is four-times larger than that occupied by glaciers (Azòcar and Brenning, 2010) and likely play an important role in the hydrological cycle (Halla et al., 2020). While melt from seasonal snow cover provides intrusive with some volcanic, volcano-sedimentary, and metamorphic rocks Paleozoic-Triassic in age (Aguilar et al., 2013;Valois et al., 2020a). Above 3000 m a.s.l. the mountain landscape has been carved by glaciers and typical geomorphology including U-shaped valleys, frontal and lateral moraines can be observed. Rock glaciers and periglacial landforms are numerous, particularly above 4000 m a.s.l.
The study site is within the semiarid Andes of Chile at the southern edge of the Arid Diagonal and Atacama Desert (Sinclair 95 and MacDonell, 2016). Specifically it is located at the headwaters of the Elqui River within the Coquimbo Region in a nature reserve called Estero Derecho (Fig. 1b). In the city of La Serena on the coast the annual precipitation is ∼ 90 mm a −1 (average from 1981-2016; Valois et al., 2020b), drastically lower than central and southern Chile. At the same time, demand from the agricultural sector, mining industry, and municipal water supply are high and water allocation has already been exhausted here (de Agua., 2016). Precipitation increases with elevation reaching ∼ 160 mm a −1 at 2900m a.s.l. in the Estero Derecho 100 valley (Valois et al., 2020b). Increased precipitation at higher altitudes allows for the formation of a seasonal snow pack that completely melts during the spring and summer seasons (Réveillet et al., 2020). Variability in precipitation at an inter-annual time scale is linked to El Niño Southern Oscillation (ENSO; Favier et al., 2009), while at a decadal time scale it is linked to the Pacific Decadal Oscillation (Núñez et al., 2013). Precipitation has decreased since 1870 by −0.52 mm a −1 at La Serena. The mean annual air temperature at a station at 3020m within Estero Derecho was 6.7°C between 2016-2020. 105 Within the nature reserve there are no debris-free glaciers, only rock glaciers and periglacial landforms. The two rock glaciers measured in this study are locally known as "El Jote" and "El Ternero" and are in the eastern part of the nature reserve (Fig. 1b) at 3700-3870 m and 4170-4510 m a.s.l., respectively. El Ternero is the largest active rock glacier within Estero Derecho, has a lobate shape and obvious flow features including ridges and furrows, a steep frontal talus slope (∼ 40 • ), and well defined lateral margins. There are a number of depressions ∼ 5 m deep on the surface and a supraglacial pond covering an area of 110 ∼ 80 m 2 . El Ternero is 1.93 km long, has a maximum width of 0.51 km, and an area of 0.60 km 2 . It is moving at a rate of ∼ 1 ma −1 and lowering by ∼ 0.2 ma −1 (based on three repeat differential GPS measurements taken in the summer of 2018-2019 and 2019-2020 between 4206 -4417 m elevation). In contrast, El Jote has poorly defined flow features and a moderately sloping frontal slope (∼ 24 • ). This landform is inactive according to repeat differential GPS measurements taken taken at five locations in the summer of 2018-2019 and 2019-2020. The lack of obvious flow features and its location within a cirque basin 115 point toward the same conclusion. El Jote is 0.86 km long, has a maximum width of 0.48 km, and covers 0.31 km 2 . Its surface is characterized by lobes as well as signs of subsidence such as depressions.
At El Ternero a stream passes adjacent to the former and eroded terminal moraine of the glacier. The waterway initiates on the mountain slope above and south of the glacier and continues down-slope, eventually feeding a high altitude wetland and the main waterway within the reserve, Estero Derecho. There is no evidence of water at the surface directly below the current 120 frontal slope of the rock glacier. However, a substantial amount of water can be heard running underneath the rock glacier within topographic depressions on the surface and it is probable that this water feeds the stream at some point below the rock glacier. At El Jote water emerges ∼ 200 m east of the main landform in a topographic low at ∼ 3740 m a.s.l. It is unclear if this water originates from the rock glacier. There is a small periglacial feature directly above that may be contributing, but no other obvious source of water visible at the surface. The waterway continues for ∼ 600 m where it disappears ∼ 100 m below 125 4 https://doi.org/10.5194/tc-2020-306 Preprint. Discussion started: 12 November 2020 c Author(s) 2020. CC BY 4.0 License. the frontal slope of the rock glacier. Water emerges in another, larger depression, along the same flow path ∼ 550 m below the front of the rock glacier and continues down-slope, contributing to an alpine wetland (i.e. bofedal) and Estero Derecho. There is vegetation adjacent to the water; in contrast there is little to no vegetation in the surrounding landscape.

Geophysics measurements 130
Surface-based geophysical methods provide information about subsurface physical properties and have been extensively used to investigate the internal structure of rock glaciers. In particular, electrical resistivity and seismic refraction tomography are common choices for the characterization of rock glacier internal structure, even though their employment in mountainous environments demands specialised techniques for sensor coupling, data acquisition and inversion routines.

135
ERT collects information about subsurface electrical resistivity (ρ) by injecting direct electric currents (DC) into the ground and measuring electric voltages at different locations. Data are obtained using a large number of resistance measurements made from spatially-distributed four-point electrode configurations (Binley and Kemna, 2005). The geometry of the current injection and potential electrode pairs are varied with typical set-ups involving many tens of electrodes and several hundred or thousand datapoints. These data are then inverted to infer the spatial distribution of electrical resistivity in the subsurface (Dahlin, 1996).

140
Electrical resistivity quantifies how strongly a material opposes the flow of electric current. In most rocks and soils, electrical current is carried by movements of ions in the pore water (electrolyte conduction), with the actual mineral matrix practically acting as an insulator (Lesmes and Friedman, 2005). Due to the high contrast in resistivity between saturated and unsaturated sediments, and the marked increase of resistivity values at the freezing point, resistivity techniques have been useful in both hydrology (de Lima, 1995;Daily et al., 1992;Valois et al., 2018b, a) and permafrost-characterization studies (Evin et al., 145 1997;Hauck et al., 2003;Langston et al., 2011). In periglacial environments, the use of ERT is particularly popular due to the contrasting electrical resistivity corresponding to lithological media, water (highly conductive) and ice (which is assumed to be an electrical insulator). In Table 1 we list the relevant values for electrical resistivity in rock glacier environments, compiled from the literature.
The main limitation for ERT is the need for the electrodes to have a good galvanic contact with the ground. Its application 150 within the surveys was therefore problematic due to the heterogeneous and rocky ground surface which resulted in extremely high contact electrical resistivity. Following Maurer and Hauck (2007) methodology, we attenuated this problem by both facilitating the injection of electric current into the ground by attaching sponges soaked in salt water to the electrodes, and in addition, increasing the measured voltage by implementing the Wenner-Schlumberger array configuration (its low geometrical factor provides larger measured voltages compared to other options). 155 5 https://doi.org/10.5194/tc-2020-306 Preprint. Discussion started: 12 November 2020 c Author(s) 2020. CC BY 4.0 License.

Refraction Seismic Tomography (RST)
RST is based on the analysis of first arrival traveltimes of critically refracted seismic waves to reconstruct seismic P-wave (i.e., compressional wave) velocity models (White, 1989;Nolet, 1987). When seismic waves impinge on velocity boundaries, they change their direction of propagation. At a critical angle that depends on the velocity contrast, head waves are created that move along the interface at the speed of the faster lower-lying layer velocity and emits refracted waves. These refracted waves 160 are measured by the receiver and the timing of their arrival (i.e., first-arrival travel times) are the main observations used in seismic refraction surveys.
Seismic velocity is the rate at which seismic waves propagate through rocks and soils and generally increase with the density of the material. In periglacial environments the large range of observed velocity values (Table 1) is favourable to the application of RST, since large velocity contrasts between the underlying materials are necessary in order to effectively interpret the data.

165
For this reason seismic refraction has been successfully used on rock glaciers since the 1970s (Barsch, 1971;Potter, 1972).
In the last two decades the method has been extensively utilized in permafrost studies (Hauck et al., 2004;Mühll et al., 2002;Draebing and Krautblatter, 2012) and to monitor hydrodynamic variation impacts on velocities (Valois et al., 2016).
One limitation of first-arrival refraction methods is that they only use a small portion of the information contained in the seismic traces and strongly depend upon the assumption that velocity increases with depth. In the case of velocity inversion 170 (i.e., the deeper medium presenting a lower P-wave velocity than the overlaying one), the refracted wave will bend towards the normal. This gives rise to the so-called "hidden layer" phenomenon (Banerjee and Gupta, 1975). Moreover, surface conditions on rock glaciers highly attenuate seismic energy and make it difficult to couple geophones and seismic sources to the ground.
During the collection of seismic data, we were able to partially improve the coupling through the use of a few geophones fastened to metal plates.We also increased the signal-to-noise-ratio by repeating and stacking the same source position multiple 175 times.

Acquisition strategy
Field data collection was conducted during the austral Summer: between the end of January and the beginning of February, 2020. The field campaign was logistically challenging in terms of both transportation of equipment and personnel to the site and difficulty moving the sensors along the profile lines because of the extremely rugged surface of the rock glaciers and the 180 altitude.
The location of sensors and sources of all the profiles were taken with a Trimble differential GPS. At both sites, we acquired the ERT surveys using Syscal Junior switch-48 (IRIS instrument, France) with 48 electrodes spaced 5 m apart and a Wenner-Schlumberger configuration. For El Jote rock glacier, the profiles length was 690 m ( Fig.1c and d), that we obtained using five sequential roll-alongs in which 50 % of the electrodes stayed in place each time and the other 50 % were displaced along the 185 profile line. In total we implemented 144 different electrode positions and obtained 2135 measurement points. For El Ternero the profile length was 575 m ( Fig.1e and f), which was obtained with four sequential roll-alongs. Here we used 120 different electrode positions and obtained 1479 measurement points. We recorded the refraction seismic surveys on both rock glaciers implementing a Geode Exploration Seismograph device (Geometric, USA) along the same lines as for the ERT profiles. The seismic source was a sledge hammer of 15 kg striking on a steel plate and we repeated each shot position five times in order 190 to improve the signal-to-noise ratio. For the profile taken on El Jote, we used 48 geophones with a spacing of 5 m and shots in between geophone positions, but spaced 10 m apart. To obtain the length of 690 m we applied five sequential roll-alongs as done for the resistivity line. In the case of El Ternero, the same spacing and configuration was used for both shots and geophones, but after the first line, the failure of one of the cables reduced the number of geophones we could use to 24. The total length of 575 m was then obtained by moving the 24-channels set-up four times and adding off-line shots to link the 195 different acquisitions. For both profiles, we manually picked the first arrival travel times on each trace resulting in 4575 picks for El Jote and 1400 for El Ternero. For both datasets the error models resulted in 1 % relative error on the ERT observations, which was obtained from the maximum of model reciprocity, and an absolute error of 0.001 seconds for the travel times, estimated from the average variability of first arrival picking. The acquisition settings are summarised in Table 2.

Joint inversion 200
Geophysical methods only provide indirect information regarding the Earth's subsurface properties, so in order to interpret them, the observations need to be inverted. With the aim of reducing the inherent ambiguity of the observations, aside from the individual inversion, we decide to jointly interpret the RST and ERT data sets by implementing two alternative approaches: structural and petrophysical joint inversion. All the inversion algorithms we use are written within pyGIMLI, an open-source library developed in python for geophysical inversion and modelling which incorporates tetrahedral/triangular unstructured 205 meshes for model discretization (Rücker et al., 2017). In each case, to quantify the models ability to explain the field observations we refer to the root-mean square error (RMSE) and the error-weighted chi-square fit, where RMSE=0 and χ 2 = 1 signify a perfect fit (Günther et al., 2006).

Structural joint inversion
Structural coupling is a flexible constraint that enforces structural similarity between the property models. The assumption 210 in this case is that the proposed models must have a similar structure since they are describing the same subsurface area, in addition to explaining the different collected observations. This similarity can be imposed either through common interfaces inversion algorithm that is suitable for 2D and 3D joint inversions on irregular meshes. The innovation in their approach is the use of a cross-gradient operator based on a correlation model acting on physical length scales, and therefore reducing the 215 dependency on a particular mesh discretization. Their cross-gradient computation is an adaptation of Levièvre and Farquharson (2013) for calculating the local spatial gradients in the model updates, where they calculate the model gradient at the location of a cell by fitting a linear trend through the model parameters of that cell and its neighbourhood. The innovation of Jordi et al.
(2019) work stands in the definition of such neighbourhood, which is given as all cells in the model that lie within a predefined distance. Their inversion scheme minimizes the following objective function: where Φ d is the combined data misfit term, λ and λ cg are the geostatistical regularization and cross gradient weights, respectively. Φ m is the geostatistical regularization (Jordi et al., 2018), which can be thought of as superposition of damping and smoothing constraints (Maurer et al., 1998). In the case of 2D problems, the main parameters to be set for this operator are the correlation lengths in the x and z directions (I x and I z ), defining the distance over which cells are considered to be correlated 225 and the angle of rotation for the correlation model (Θ) with respect to the Cartesian plane, so that the rotated coordinate axes can be set parallel to the geologically dominant directions. The last term, Φ cg , is the cross-gradient operator and it is defined on the same scale as the regularization term.
In order to apply the structural joint inversion scheme, we first verify that the results from individually inverted data-sets were not completely different, so that there is a structural similarity (Linde et al., 2008).
The treatment of the rock volumetric fraction as a single phase is a justified simplification in rock glacier environment, where the amount of soil is negligible compared to the hard rock.
The volumetric fractions in Eq.
(2) are related to the seismic slowness (s), reciprocal of the P-wave propagation velocity (v), 240 through the time averaging equation (Timur, 1968): and to the electrical resistivity through a modification of Archie's second law (Archie, 1942): where the porosity is expressed in terms of rock content (φ = 1 − f r ). The assumptions within this 4PM model are that the 245 medium is isotropic, at high effective pressure and has a single homogeneous mineralogy (validity of Eq. 3), and that the electric current flow is dominated by electrolyte conduction (validity of Eq. 4).
The petrophysical joint inversion scheme minimizes the following objective function (Wagner et al., 2019;Mollaret et al., 2020): Similarly to Eq. (1), the first term refers to the combined data misfit, while Phi m represents a smoothness regularization term build through four first-order roughness operators to promote smoothness in the distribution of each constituent of the four-phase system. The last term is an additional regularization term which constrain the volume conservation (Eq. 2). The two weights λ and λ p are responsible for scaling the influences of the two regularization terms, where λ is chosen to fit the data within the error bound and λ p is chosen large enough to prohibit no-physical solutions.
Within this framework, the RST and ERT observations are used to infer the volumetric fractions of water, ice, air and rock for each model cell, while the spatial distribution of electrical resistivity and P-waves velocities are obtained through Eq. (3) and Eq. (4), where the petrophysical parameters and constituent velocities are assumed to be spatially constant. We chose the values for the inversion of the field observations based on the literature (Hauck and Kneisel, 2008;Maurer and Hauck, 2007;Hauck et al., 2011;Wagner et al., 2019), which are listed in Table 3. A last important parameter to consider in this scheme is 260 porosity initial value and range. Wagner et al. (2019), already stressed the importance of a good porosity estimation in order to avoid ambiguity between ice and rock content and in a recent study, Mollaret et al. (2020) analyses the influence of the porosity constraint in the petrophysical joint inversion results. Following the approach of this last paper, and accordingly to the previous knowledge from the field site we tested different initial porosity values and ranges (φ min -φ max ) for both rock glaciers. The choice was made selecting the less constraining intervals which allowed results consistent with the geology of the two sites.

El Jote
Within the profiles taken on El Jote, the electrode positions coincide with the geophone locations. To interpret the observations obtained though seismic and resistivity tomography, we implemented individual dataset inversion and both petrophysical and structural joint inversion schemes as described in the Methods section. In all cases, we set a homogeneous resistivity starting 270 model, with a value equal to the median of the apparent resistivities (ρa median = 4561 Ω m) and a gradient model for the seismic velocity, starting with 300 m s −1 at the top of the domain and gradually increasing to 5000 m s −1 at the base. The regularization weight was set to λ = 10 for all inversion schemes, chosen according to L-curve analysis (Hansen, 2001). As introduced in the methodology session, we firstly invert individually the data-sets in order to ensure structural similarity. The model results for El Jote are given in Fig. 2(a) and (b).

275
For the structural joint inversion scheme, we set the vertical and horizontal correlation lengths for the geostatistical regularization to I x = 40 m and I z = 15 m, respectively, and the angle of rotation to Θ = −18.8 • (equal to the dominant profile inclination). In addition, we set the cross gradient weight to λ cg = 1000 after L-curve analysis of the joint χ 2 values. We obtained an overall data fit of χ 2 = 1.49 after 15 iterations with RMSE for the ERT and traveltime data of 1.62 % and 1.55 %, respectively. Figure  In Fig. 5 we compared the relative residuals of the three inversion schemes computed as: with d OBS standing for the field observations (ERT or RST) and d M OD for the response of the inversion model results.

El Ternero
Collection of the profiles on El Ternero were logistically more complicated than on El Jote, due to higher altitudes and especially to the extremely heterogeneous, rocky surface. Even though resistivity and seismic lines were identical, some variation 310 can exist between geophone position (on a rock in case of flat geophone) and electrode position (between rocks). The position differences are much lower than the spacing between two sensors, so that such positioning error should not have a significant impact on the results. In order to implement the joint inversion schemes we therefore needed to adapt the sensor locations to lie within the same line. Because RST is more sensitive to the local geometry than ERT, we decided to adapt the electrode locations to the geophones line while preserving the relative distances of the electrodes. Overall, the data quality for this rock 315 glacier is much lower than for El Jote, which is impacted by the significant decrease in the number of data points for both ERT (almost 1.5 times less than for El Jote) and RST (more than 3 times less than for El Jote, see Table 2). As for El Jote, we initiated the inversion schemes with a homogeneous resistivity, with values equal to the median of the apparent resistivities (ρa median = 36054 Ω m) and the same gradient model for the seismic velocity (i.e., 300 m s −1 at the top and 5000 m s −1 at the bottom of the domain). We applied the same regularization and cross-gradient weights as for the other site, only changing the angle of rotation for the regularization operator to the dominant profile inclination for El Ternero: Θ = −10 • . As for the previous case, we firstly individually invert the observation to asses structural similarity (see Fig. 2c and d).
From the structural joint inversion, we obtain an overall data fit of χ 2 = 1.51 after 15 iterations with RMSE for the ERT and traveltime data of 1.90 % and 0.91 %, respectively. In Fig. 6  For the petrophysical inversion scheme, the initial porosity was set to an homogeneous 60 % and inverted within a range going from 10 % to 90 %, heterogeneously within the model. After 15 iterations we obtained an overall data fit corresponding 335 to χ 2 = 1.65 with RMSE for the ERT and traveltime data of 1.78 % and 1.30 %, respectively. In Fig. 7 the joint inversion results obtained through petrophysical coupling are shown, which help to clarify the information gained through the structural joint inversion. Indeed, they confirm the presence of a thin top layer with moderately high fraction of air (up to 28 %, see 17 %, respectively (Fig. 7c). As in the previous case, the transformed velocity and resistivity models ( Fig. 7a and b) present differences if compared to the structural joint inversion, with overall lower velocity values and higher resistivities values and contrasts. In Fig. 8 we compared the relative residuals of the three inversion schemes for El Ternero glacier (Eq. 6). For both rock glaciers the structural coupling is able to slightly increase the data fit for RST observations, while the petrophysical coupling presents an improved ERT data fit (Table 4). In both cases, the RMSE obtained in the case of joint inversion increases only slightly if compared to the inversion of individual datasets, proving that it is possible to decrease the ambiguity of the inversion results without significantly losing accuracy for the model responses in describing the observations.

375
For El Jote, the results show a top layer (laterally variable between 10 to 50 m thick) of unconsolidated rock with air-filled pore space, especially from 300 m to the end of the profile line. This overlays a layer where the pore space decreases and appears saturated with water for the majority of the line, apart from near 550 m, where the fraction of ice slightly increases to 3 % (Fig.4d). The increased velocity and resistivity at 550 m could also be interpreted as the presence of intact rock, as opposed to ice. Indeed, when the porosity of the subsurface is unknown, the petrophysical joint inversion scheme does not easily 380 differentiate between ice and rock content (Wagner et al., 2019). In order to gain information about porosity, we unsuccessfully attempted to drill a core sample, but the drill broke at very shallow depths. Nevertheless, for both interpretation outcomes, the overall results are a reflection of the permafrost degradation on El Jote. It is likely that within this inactive rock glacier, the ice has melted, leaving behind large voids filled with air (top layer) or water (deeper layer). For both the inversion results it seems that the bedrock is deeper than 100 m for almost the entire profile length, but the strong increase in velocity (Fig. 3a) and rock 385 content (Fig. 4f) at the bottom of the domain and at profile lengths of 150 m to 250 m and at approximately 550 m may be interpreted as the top of the bedrock. In addition, the high water content in the bottom layer (more than 20 %) suggests the presence of an aquifer between the bedrock and the surface of the inactive rock glacier.

El Ternero (active rock glacier)
The inversion model results for El Ternero are slightly shallower than the ones obtained on El Jote. This is due to the extremely

Towards a diagnostic model representation for the ice or water presence in rock glaciers
Both joint inversion schemes give satisfactory results which are overall consistent with each other for El Jote, but leave some ambiguity in case of El Ternero. Nevertheless these interpretation methods are not always applicable as the petrophysical joint inversion may be limited by the lack of proper petrophysical models (or parameters), and the structural coupling may 405 not be suitable due to the inherent properties of the different data sets ( e.g., resolution, signal-to-noise ratio and sensitivity to lithological structures). When model coupling is not possible, the comparison of velocity and resistivity model inversion results can still deliver plenty of information about the rock-glacier's internal structure. For example, in Fig. 9 we show scatter plots of resistivity vs velocity for (a) El Jote and (b) El Ternero individual (blue asterisks) and structural joint inversion results (orange asterisks). In the first case, we can see for both inversion approach a bulk of points characterized by velocity lower 410 than 1000 m s −1 and relatively high resistivities (ρ ∼ 10 4 Ω m) which correlates with the large air content modelled at the top of the investigated domain by petrophysical coupling (Fig. 4e). The overall trend is represented with a solid line (red for individual and black for structural joint inversion results), which has a negative coefficient of slope which is larger for the structural joint inversion (m=-0.28) than for the individual inversion (m=-0.18). This means a decrease in resistivity (reduction in air or increase in water content) while velocity increases (reduction of the pore space or consolidation with ice presence 415 or healthier rocks). This trend is coherent with the hypothesis of a subglacial aquifer, already deducted by the water content 13 https://doi.org/10.5194/tc-2020-306 Preprint. Discussion started: 12 November 2020 c Author(s) 2020. CC BY 4.0 License. obtained by petrophysical joint inversion (Fig. 4c). For El Ternero, the scatter plot shows a rather different behavior: a large density of points presenting high resistivity (ρ > 10 4 Ω m) and velocities between 1000 and 3500 m s −1 , followed by a general trend of increasing resistivity and velocity. Also in this case, the slope coefficient in case of structural joint inversion (m=0.43) is larger than the one from individual inversion (m=0.27). Such a trend is consistent with the ice rich layer depicted in Fig.   420 7(d), with higher resistivities and velocities. Moreover, for both scatter plots the results from individual inversion are more disperse compare to the one from structural joint interpretation. The rather different appearance of the two scatter plots can be used as an indicator of the distinct nature of the two rock glaciers: overall, the inactive rock glacier is characterized by lower resistivities and a large range of velocities representing either air (v < 1000 m s −1 ) or water filled pore space, while the active rock glacier is indicated by higher resistivity and velocity values, reflecting the ice rich layer. Moreover, in the case of the active 425 rock glacier the slope coefficient is positive, while it is negative in the inactive case. This type of model representation is a valuable tool for interpreting the inversion model results and can be implemented as well in case it is not possible to implement neither structural nor petrophysical coupling.

Hydrogeological role
El Ternero and El Jote represent two end-members of rock glacier types. El Ternero is an active rock glacier containing 430 significant amounts of ice according to our geophysical analysis (Fig. 7d), while El Jote is likely an inactive rock glacier whose ice has largely if not completely melted to form an aquifer ( Fig. 4c and d). Each have a distinct and important hydrological role. El Ternero has the capacity to function as a long-term aquifer given that most of the water it contains is in the form of ice which is insulated from the environment by debris cover (∼ 5 m thick). The insulating effect of the debris cover has been shown to slow the rate of melt (Jones et al., 2018;Bonnaventure and Lamoureux, 2013) making rock glaciers more resilient to 435 climate change compared to debris-free glaciers. Field observations and results from the joint inversion modeling suggest that some of the ice contained within El Ternero is melting. There are various oval-shaped topographic depressions ∼ 5 m deep on the surface of the rock glacier including a supraglacial lake that is ∼80 m 2 . We interpret these depressions as thermokarst degradation features. In some of these depressions a substantial amount of water can be heard flowing which may be derived from melting of ice within the rock glacier. The joint inversion results show that there may be water accumulating at 450 m 440 along the profile (Fig. 7c).
El Jote has an important hydrological role as an aquifer and in terms of its influence on water flow. Assuming its impact on hydrology is similar to other rock glaciers for which measurements exist in the semiarid Andes and elsewhere, it is likely to act as a reservoir during high melt and/or precipitation periods then release the water downstream at a rate slow enough to contribute significantly to streamflow when it is needed most at the end of summer (Pourrier et al., 2014;Harrington et al., 445 2018). It is possible that the surface streams down-glacier are the surfacing of water from the rock glacier aquifer, although this hydrological link has not been measured in the field. This role in water storage and delayed runoff is critical for the Coquimbo Region where annual precipitation is very low (<100 mm a −1 at the coast), precipitation events are variable, and there is little to no precipitation during the warmest months of the year (Garreaud, 2009;Valois et al., 2020a).
In this study, we presented the first comparison of geophysical signatures of one active and one inactive rock glaciers using independent, structurally coupled and petrophysical coupled inversion routines. We found that each routine agrees with much higher velocities and resistivities for the active rock glacier, which are interpreted as a much higher content of ice in the petrophysical model. Scatter plots of velocity versus resistivity point out as well a clear signature difference for both rock glacier, which could be interpreted as the difference between a relict glacier with an aquifer and an active glacier containing 455 ice-rich layers.
The analysis of the inversion results confirms the effectiveness of geophysics for rock glacier (active or inactive) characterization but at the same time proves the necessity for complementary geophysical measurements: only velocity or resistivity information does not allow the differentiation between the distinct subsurface components (i.e., air, water, rock and ice). The

460
proving that it is possible to decrease the ambiguity of the inversion results without significantly losing accuracy for the model responses in describing the observations.Both joint inversion schemes could be used on their own to interpret the observation, nevertheless the comparison between the results of both schemes helped us to reduce some of the model results ambiguity.
Through the joint interpretation of ERT and RST surveys for El Jote we were able to characterize its inner structure, detect the top of the bedrock in part of the model domain and identify a potential aquifer, while in case of El Ternero we could 465 clearly determine the active layer and the top of the permafrost, together with the sign of its partial melting at the bottom of the sensed domain. Nevertheless,there is ambiguity in the interpretation between ice and rock matrix (especially for El Jote inversion results), which could be improved adding information about subsurface porosity or by the incorporation of additional freeze-thaw sensitive data sets such as complex electrical resistivity measurements (Wagner et al., 2019). In addition, to increase the investigated depth it would be necessary to improve the seismic data quality, which could be (time-expensively) 470 done by fastening the geophones to the surface by drilling small holes in the rock.
tremendous help in the field. We also wish to thank Claudio Jordi for sharing his code and Florian Wagner and his research group for making it available to the public. We thank Ignacio Castro Cancino for his contribution to the rock glacier descriptions and for providing the MAAT data for the station in Estero Derecho. This work was supported by CONICYT-Programa Regional-Fortalecimiento (R16A10003) and FIC-R  29 https://doi.org/10.5194/tc-2020-306 Preprint. Discussion started: 12 November 2020 c Author(s) 2020. CC BY 4.0 License.  Kneisel, 2008 andHauck, 2007).