Introduction
Accurate ice surface displacement measurements are important inputs to the understanding of ice dynamics and changes in mass balance of glaciers. Ice flux is largely dependent on glacier velocity, with accelerated or decreased motion indicating an alteration of ice mass equilibrium (Reference Berthier and VincentBerthier and Vincent, 2012; Reference Nuth, Schuler, Kohler, Altema and HagenNuth and others, 2012) or ice dynamics (Reference Sakakibara and SugiyamaSakakibara and Sugiyama, 2014). Reference Heid and KääbHeid and Kääb (2012) observe that the decreasing surface velocity of land-terminating glaciers over several years is correlated with negative mass balance reflecting lower mass transport from the accumulation to the ablation areas. In contrast, calving glacier velocities can increase despite periods of negative mass balance (Reference Rivera, Bown, Carrión and ZentenoRivera and others, 2012; Reference Sakakibara, Sugiyama, Sawagaki, Marinsek and SkvarcaSakakibara and others, 2013; Reference Sakakibara and SugiyamaSakakibara and Sugiyama, 2014). Ice surface velocity measured over periods of a few years, in combination with the principle of mass continuity, can also be used to derive ice thickness estimates (Reference McNabbMcNabb and others, 2012; Reference Farinotti, Corr and GudmundssonFarinotti and others, 2013). Measuring glacier surface velocities at shorter timescales has proved crucial to understanding the link between glacier hydrology and basal motion (Reference Iken and BindschadlerIken and Bindschadler, 1986; Reference Hooke, Calla, Holmlund, Nilsson and StroevenHooke and others, 1989; Reference Iken and TrufferIken and Truffer, 1997; Reference Sugiyama and GudmundssonSugiyama and Gudmundsson, 2004; Reference SugiyamaSugiyama and others, 2011). Among other applications, surface velocity at daily to hourly timescales can be used to improve understanding of the forcing and mechanisms behind surge events and/or floods in glacier-fed rivers (Reference Bartholomaus, Anderson and AndersonBartholomaus and others, 2008).
Differential global navigation satellite system (GNSS) ground surveys are the most accurate method to obtain point-wise glacier velocity measurements as they can resolve sub-seasonal to hourly velocity fluctuations. However, in remote areas such as the North Patagonian Andes it is difficult to perform regular ground-based surveys over glaciers. In addition, point GNSS measurements cannot (unless we have many more measurements available) depict the spatial pattern of ice surface velocities over the entire glacier surface and it is therefore necessary to use indirect alternative techniques. With the development of automatic and semi-automatic remote-sensing techniques such as image-to-image cross-correlation (Reference Bindschadler and ScambosBindschadler and Scambos, 1991; Reference Leprince, Barbot, Ayoub and AvouacLeprince and others, 2007), it is now possible to retrieve glacier surface velocities by matching multitemporal satellite images and this technique has become a well-established methodology in glaciology (Reference Heid and KääbHeid and Kääb, 2012).
Most studies of glacier velocity changes in the Patagonian Andes have focused on the large outlet glaciers of the South and North Patagonian Icefields, and the Cordillera Darwin Icefield in Tierra del Fuego (Reference Naruse, Fukami and AniyaNaruse and others, 1992; Reference Michel and RignotMichel and Rignot, 1999; Reference Skvarca, Raup and De AngelisSkvarca and others, 2003; Reference Melkonian, Willis, Pritchard, Rivera, Bown and BernsteinMelkonian and others, 2013; Reference Sakakibara, Sugiyama, Sawagaki, Marinsek and SkvarcaSakakibara and others, 2013). The only surface velocity data available for glaciers in the North Patagonian Andes were reported by Reference Rivera, Bown, Carrión and ZentenoRivera and others (2012) for Los Amarillos glacier on the flank of Michimahuida volcano (∼42.8° S, 72.4° W). Those authors studied the dynamic changes of this glacier, between 2007 and 2011, noting an increase in the surface velocity of the lower reaches associated with an advance of the glacier in 2009–11. The behavior of this glacier contrasts with the general recession of glaciers in the North Patagonian Andes, and was tentatively related to the volcanic activity in the area (Reference Rivera, Bown, Carrión and ZentenoRivera and others, 2012).
In this paper, we apply cross-correlation to high-resolution images acquired by the Pléiades satellite (Reference Gleyzes, Perret and KubikGleyzes and others, 2012) to generate high-resolution sub-seasonal surface velocity maps of glaciers of Monte Tronador in the North Patagonian Andes of Argentina and Chile. These maps have been used to describe the spatial flow patterns and magnitude of surface velocities on the glaciers and to analyze the spatial and sub-seasonal changes in relation to potential inputs of meltwater and/or rain to the subglacial system. To our knowledge, although Pléiades images have previously been used to map surface displacements of landslides (Reference Stumpf, Malet, Allemand and UlrichStumpf and others, 2014; Reference Lacroix, Berthier and MaquerhuaLacroix and others, 2015), this is their first use to derive surface velocities on glaciers.
Study Area
Monte Tronador (3475 m a.s.l.; 41.15° S, 71.88° W) is an extinct stratovolcano located in the North Patagonian Andes along the Argentina–Chile border (Fig. 1). The upper slopes of Monte Tronador host the accumulation zones of eight major valley glaciers, namely Frías, Alerce, Castaño Overa and Manso/Ventisquero Negro in Argentina, and Verde, Blanco, Norte and Casa Pangue, as well as three unnamed glaciers, in Chile. Casa Pangue, Verde and Manso/Ventisquero Negro have regenerated debris-covered tongues below steep icefalls/cliff walls at ∼1400–1500 m a.s.l. All other glaciers are largely debris-free and do not descend below the 1500–1400 m a.s.l. bedrock step.
Fig. 1. False-color pansharpened Pléiades image of Monte Tronador from 7 March 2012 (RGB bands 3, 2 and 4) showing the individual glaciers (delineated in blue) and other locations mentioned in the text. Red dots locate the ablation stakes where velocities were measured in the field. Light blue indicates the flow profiles discussed in Figure 9. Rectangles identify Figure 2a and c. The longitude and latitude of the summit of Monte Tronador are shown in white for reference. The right and upper grid axes show UTM 19 S coordinates, and the left and bottom grid axes show the scale in meters. The same projection and scale display is used for other maps. © CNES 2012, Distribution Airbus D&S.
The seasonal variation in precipitation and weather in the North Patagonian Andes is largely modulated by the location and intensity of the Southern Hemisphere westerlies (Reference GarreaudGarreaud, 2009). The equilibrium-line altitude (ELA) at Monte Tronador is ∼2000 m a.s.l. (Reference Carrasco, Casassa and QuintanaCarrasco and others, 2005, Reference Carrasco, Osorio and Casassa2008; Reference Condom, Coudrain, Sicart and ThéryCondom and others, 2007) and numerical models estimate that these glaciers have a high mass turnover rate with maximum accumulation and ablation of more than +5 and −10 m w.e. a−1, respectively (Reference Leclercq, Pitte, Giesen, Masiokas and OerlemansLeclercq and others, 2012). Our preliminary mass-balance measurements, conducted in 2013–14, confirm these values.
All glaciers in this area were substantially larger, developing clear moraines and trimlines during the Little Ice Age (LIA) between about 1650 and 1850 (Reference Masiokas, Rivera, Espizua, Villalba, Delgado and AravenaMasiokas and others, 2009, Reference Masiokas, Luckman, Villalba, Ripalta and Rabassa2010). Presently, however, most of these glaciers show a clear retreating and/or thinning pattern (Reference BownBown, 2004; Reference Masiokas, Rivera, Espizua, Villalba, Delgado and AravenaMasiokas and others, 2009; Reference Paul and MölgPaul and Mölg, 2014). Casa Pangue and Río Manso glaciers have experienced significant thinning (Reference Bown and RiveraBown and Rivera, 2007; Reference Masiokas, Luckman, Villalba, Ripalta and RabassaMasiokas and others, 2010), and all glaciers show a long-term frontal recession since the beginning of the 20th century except for Verde glacier, which remains in contact with LIA moraines. Our observations indicate that debris-covered glaciers such as Verde, Casa Pangue and Río Manso have experienced lower retreat rates than the debris-free glaciers. Reference Bown and RiveraBown and Rivera (2007) and Reference Masiokas, Villalba, Luckman, Lascano, Delgado and StepanekMasiokas and others (2008) suggest that a regional warming trend, together with a concurrent decrease in precipitation, partly explains the recent regional retreat observed in North Patagonian glaciers.
Reference Villalba, Leiva, Rubulls, Suarez and LenzanoVillalba and others (1990) and Reference Masiokas, Luckman, Villalba, Ripalta and RabassaMasiokas and others (2010) dated the LIA maximum moraines of Frías (clean ice tongue) and Manso/Ventisquero Negro (debris-covered tongue) glaciers to ∼1600 and the mid-1800s, respectively. Reference Masiokas, Luckman, Villalba, Ripalta and RabassaMasiokas and others (2010) tentatively concluded that the asynchronous responses of these glaciers were related to differences in valley and glacier morphology and dynamics. However, this hypothesis was difficult to validate due to the lack of detailed glaciological and climatological data. During the 1970s Reference Rabassa, Rubulis and SuarezRabassa and others (1978) reported a small but rapid advance of all glaciers on the Argentinian side of Monte Tronador. The advance of Frías glacier was so fast that they hypothesized the possible occurrence of a surging event.
Materials and Methods
The analyses and results presented in this study are based on five images acquired by the Pléiades-1A satellite between 7 March and 6 June 2012 (Table 1). These images are delivered at a ground sampling distance (pixel size) of 0.5 m for the panchromatic channel (480–830 nm) and 2 m for the multispectral blue, green, red and near-infrared channels (http://smsc.cnes.fr/PLEIADES/). The actual resolution of the sensor, however, is slightly coarser (0.7 and 2.8 m) and the images are therefore oversampled by the ground segment (Reference Gleyzes, Perret and KubikGleyzes and others, 2012). Pléiades panchromatic band images are coded on 12 bits (i.e. 4096 digital numbers), which gives much better image contrast than the classical 8-bit encoding. This is particularly advantageous for the identification of features over flat and textureless regions such as the snow-covered accumulation areas of glaciers (Reference BerthierBerthier and others, 2014).
Table 1. List of satellite images used in this paper. Pléiades-1A satellite (DS_) covers the period between 7 March and 6 June 2012. Spatial resolution values are for the original data and (in parentheses) the orthoimages generated in this study
The following sections describe the procedures undertaken to generate ice surface velocity maps from the raw (level 1A) satellite images. A description of the degree-day model used in the discussion section and its calibration is summarized in the Appendix.
Orthorectification of images
The orthorectification of Pléiades images was conducted using the Orthoengine module of the software PCI Geomatica 2013. The original orientation of each image was obtained from the rational polynomial coefficients available in the ancillary data. Without ground control points (GCPs), the horizontal location accuracy of the images was estimated at 8.5 m (CE90, circular error at a confidence level of 90%; Reference Lebègue, Butler, Xiong and GuLebègue and others, 2013). The initial orientation was then refined using five GCPs and between 20 and 30 tie points collected in the overlapping region of the five images. The GCPs were measured in the field with a Trimble GeoExplorer 6000 Series GNSS receiver in static mode. Each GCP was measured for at least 45 min, and differential correction was performed using the Osorno permanent GPS station (84 ma.s.l.; 40.59°S; 73.10°W) located ∼120 km west of Monte Tronador. The horizontal and vertical precision of the differentially corrected coordinates of the GCPs are 0.1 m and 0.2 m, respectively. Table 1 lists the accuracy of the GCPs and tie points for each image. Finally, 1 m spatial resolution, panchromatic orthoimages were generated using a DEM derived from a triplet of Pléiades stereo images acquired on 21 April 2012. The DEM was generated following the procedure described in Reference BerthierBerthier and others (2014), and has a spatial resolution of 2 m and an accuracy (root-mean-square error (RMSE)) of 0.5 m and 1.1 m in the horizontal and vertical directions, respectively.
Cross-correlation of orthoimages
Horizontal ground displacements were retrieved from the sub-pixel correlation of successive Pléiades orthoimages using the COSI-Corr package, a plug-in for the Environment for Visualizing Images (ENVI) software (Reference Leprince, Barbot, Ayoub and AvouacLeprince and others, 2007). Once orthorectified and co-registered, the correlation between two images results in two offset fields, representing the east–west and north–south horizontal ground displacements, and a signal-to-noise ratio (SNR) field reflecting the quality of the correlation for each pixel in the images. To allow the measurement of large displacements without losing resolution on the displacement fields, COSI-Corr multiscale correlation analysis was performed using an initial search window size of 256 down to 8 pixels with steps of 16 pixels. The initial displacement is first estimated based on search windows, and is then refined in the subpixel domain. An important condition is that the range of window size has to be large enough to capture the expected higher displacement. In our case, with a maximum window size of 256 pixels, we could measure a maximum surface velocity between 9 and 14 m d−1, which is seven to nine times larger than the measured maximum surface velocity (1.2 m d−1). Finally, the surface velocity maps are sampled every 16 m for Pléiades. Table 2 indicates the pairs of Pléiades images cross-correlated in this study.
Table 2. Image pairs used to obtain surface velocity measurements for the glaciers of Monte Tronador. ID indicates the informal name for each pair of images used in the text
Generation of surface velocity maps
We applied a series of filters to the output of the cross-correlation analyses to exclude poorly correlated pixels and those showing unrealistic displacement values and/or anomalous flow directions. We initially excluded all pixels with SNR values lower than 0.95 to avoid miscorrelation and subsequently applied a magnitude threshold filter. We examined the regions of coherent high velocities to determine the maximum realistic velocities over glaciers in each period. The maximum coherent surface speed was ∼1.2 m d−1 in all periods. Values above this limit were discarded on and off glaciers before the accuracy analysis. Secondly, we applied a flow direction filter over glaciated regions (Reference Berthier, Raup and ScambosBerthier and others, 2003). Measurements that deviated by more than ±60° from the aspect of the slope (obtained from the Pléiades DEM) were excluded. This sequential application of filters allowed for efficient identification and removal of most anomalous values due to poor correlation, shadows, avalanches or serac falls, the upward tilt of seracs, and changes in the snow cover (Fig. 2a–d).
Fig. 2. The influence of the successive filters applied to the velocity field. The panels compare the image pairs used to derive the March velocity field. Date format is dd/mm/yyyy. (a, b) The serac and avalanche falls of Manso glacier. Note how the cells located in the avalanche cone pass the SNR filter but not the directional filter, and how most of the cells representing the debris-covered tongue pass all of the filters. (c, d) The bedrock outcrop that separates Frías and Manso glaciers. Most of the values associated with the change in the snow cover of the bedrock outcrop do not pass the directional or magnitude threshold filters, but the glacier-covered areas show a coherent pattern of flow.
Finally, glacier surface velocity maps for each period (March, April, May and June; Table 2) were obtained after applying a spline interpolation to fill data gaps generated by the data exclusion in previous steps. A spline interpolation results in a smooth surface without modifying the original values of pixels adjacent to the voids when they are small, but can results in spurious undulations over large voids. Thus, we only interpolated data gaps smaller than 0.15 km2 and left larger gaps (only on the June surface velocity map) unfilled. Finally, we estimated the mean annual surface velocity for glaciers of Monte Tronador for 2012, from the mean displacement measured between March, April and May.
Accuracy of displacement measurements over ice-free areas
The accuracy of the ice surface velocity measurements was first assessed by a null test over motionless, ice-free areas (Reference BerthierBerthier and others, 2005). To separate glaciated from ice-free areas, we manually digitized a glacier mask from the 27 March 2012 panchromatic Pléiades orthoimage, the one with least seasonal snow cover.
Several studies have used the mean, standard deviation and/or RMSE over ice-free areas as an indication of displacement errors based on the assumption that errors of the displacement measurements are normally distributed (Reference BerthierBerthier and others, 2005; Reference Scherler, Leprince and StreckerScherler and others, 2008; Reference Scherler and StreckerScherler and Strecker, 2012). Reference Höhle and HöhleHöhle and Höhle (2009) present a method to assess whether the error (in our case the displacements on the ice-free terrain, hereafter referred to as the ‘residuals’) is normally distributed. If it is not, they propose an alternative robust metric statistic to quantify the residual populations. Their method is primarily intended to be used with digital elevation models (DEMs) but it can be used to measure the accuracy of other residual populations with heavy tails. Those authors propose first to use a Q–Q plot (comparing data with a statistical population) to intuitively ‘see’ whether the distribution is normal and, if not, to use quantiles as a measurement of error because they are less sensitive to large outliers. Reference Höhle and HöhleHöhle and Höhle (2009) propose to use the median (Q 0.5) and the quantiles 0.95 and 0.683 of the absolute values of the residuals to quantify the error in the measurements.
We used the normal probability plot to assess the distribution of the residuals. In this technique, the plot of residual data against a theoretical normal distribution should approximate a straight line. The results indicate that the distribution of the residuals was not normal, with heavy tails (Fig. 3) for all periods under investigation. Thus, we used the median Q(0.5), and the quantiles Q(0.683), Q(0.9) and Q(0.95) of absolute residuals to quantify the error in the displacement measurements (Table 3). The distribution of residuals for each period (Fig. 4) shows that the value of Q(0.683) excludes a high percentage of data as outliers, while Q(0.95) includes as data residuals that are related to the displacement in areas in shadow (Fig. 5). We found that Q(0.9) was a good compromise to separate most outliers (especially those associated with displacements in shadows) from valid data obtained in illuminated sectors (Figs 4 and 5).
Fig. 3. Normal probability plot for the pair of images for April that shows that the error distribution is not normal. The absolute displacement values (blue dots) considerably deviate from a linear relationship (red dashed line).
Fig. 4. Absolute displacement measurements on the ice-free terrain for each Pléiades image pair, showing the values of the different quantile values that could be used to measure the error (see Table 3). Date format is dd/mm (2012).
Fig. 5. (a) Relationship between the total number of data points and outliers based on Q(0.9) and aspect. Note a concentration of outliers in the southeast to southwest azimuth. (b) Same as (a) but with slope. In this case, outliers are concentrated around a slope of 40°. (c–f) A perspective view of Monte Tronador and the distribution of outliers (colored) using the different quantiles as thresholds. Note how the Q(0.683) includes most of the data in forested areas as outliers, whereas the outliers of Q(0.90) are restricted mainly to cast-shadow areas. Q(0.95) does not include as outliers some of the pixels that are in shadows.
Table 3. Error measurements (Q(0.90)) for each period assessed in this study. The last column shows the error value using Q(0.9) expressed as surface velocity
It is important to note that, because of the high spatial resolution of the Pléiades images, the displacement over most of the non-glaciated areas (Table 3) is likely affected by change in the cast shadow of the tree canopies (e.g. Reference BerthierBerthier and others, 2005; Reference Stumpf, Malet, Allemand and UlrichStumpf and others, 2014). This effect depends on the date and time of the scenes and cannot be corrected. However, this effect has no influence on glaciers, suggesting our error estimates are likely conservative. We selected areas without trees or shadows to calculate Q(0.9) and test this hypothesis.
Accuracy of displacement over glaciers
The lack of ground-truth surface velocity measurements during the same time span covered by the satellite images hampers a simple evaluation of the remotely sensed glacier displacements obtained in this paper. Thus, we assessed the overall quality of these measurements using different approaches. We constructed streamlines from the displacement fields generated and compared their direction with flow features on the glaciers in the orthoimages. We also compared the obtained magnitude and direction of displacements with precise GPS measurements collected at four ablation stakes between March and September 2013 (see Fig. 1 for location). After differential correction, the uncertainty for the displacement measurements at these stakes was 0.14 m.
Results
Residuals over ice-free areas
The surface velocity errors are between 0.03 m d−1 (11 m a−1) during March and April, 0.06 m d−1 (22 m a−1) in May and 0.28 m d−1 (102 m a−1) in June (Fig. 4; Table 3). We believe that this increase of one order of magnitude in the error in June is related to the poor cross-correlation between images due to the increase in snow cover from 12 May to 7 June (Fig. 4c; Table 4). We found that most outliers are located in areas facing south-southeast–southwest (150° to 270°), opposite to the sun azimuth angle for the scenes (between 30° and 40°; Table 1). This indicates that outliers tend to be spatially correlated with shadow zones (Fig. 5). We also found that error is greatly reduced (up to 50% lower in all periods except for the June pair; Table 3) in areas without trees or shadows. This confirms our hypothesis that our error estimates are conservative and affected by the displacement of the cast shadows of the tree canopy in ice-free areas. Nevertheless, these areas without trees or shadows represent <0.5% of the ice-free terrain and had to be selected manually, which could introduce other biases. So, conservatively, we keep the Q(0.9) for the whole ice-free terrain as our error estimation (Table 3).
Table 4. Morphometric parameters and mean annual surface velocity for the glaciers of Monte Tronador
Surface velocity maps
Figure 6 shows the surface velocity maps obtained from each Pléiades pair before gap filling, and Figure 7 shows the mean annual surface velocity map for glaciers of Monte Tronador. Table 4 presents the mean surface velocities for individual glaciers of Monte Tronador and for the different periods investigated in this study. These surface velocity maps confirm the radial flow pattern from the peak on Monte Tronador. All periods have a similar spatial flow pattern over the glaciers of Monte Tronador, justifying combining them in a synthetic velocity map. The highest velocities occur in the steep icefalls of Verde, Casa Pangue and Manso glaciers and exceed 350 m a−1. The gently sloping Alerce and Blanco Chico glaciers have the lowest velocities (mean of 16–60 m a−1). Several fine-scale glacier flow characteristics are well represented by the velocity maps (e.g. the increase in velocity in the funnel zones of Norte and Blanco Chico glaciers, the decrease in velocity in the lower reaches of the debris-covered tongues of Verde and Casa Pangue, and even the change of surface velocity associated with minor slope change in the accumulation areas of Castaño Overa and Alerce glaciers). Also, the maps reveal a clear acceleration in the lower parts of Frías and Manso/Ventisquero Negro glaciers. They also show that ice flow is parallel to some of the ridges (Fig. 7).
Fig. 6. Surface velocity maps for each study period. Date format is dd/mm (2012).
Fig. 7. Estimated annual surface velocity of glaciers of Monte Tronador for 2012. The data for total surface displacement measured between 7 March and 12 May 2012 (66 days) were extrapolated to the rest of the year. White arrows show the mean surface direction aggregated for 15 neighboring cells, which cover ∼0.05 km2. The inset shows the ice divide between Frías and Casa Pangue glaciers, with the flow direction arrows representing the mean direction of five neighboring cells (0.006 km2).
Comparison with displacement of ablation stakes
Figure 8 shows the comparison between velocities measured at ablation stakes and those calculated from Pléiades velocity maps. Ablation stake S.4 is located close to the ice divide between Castaño Overa and Alerce glaciers, where the ice motion is very low without sub-seasonal variation. During both periods, the GNSS and the Pléiades velocity measurements are very similar. At ablation stake S.9, located in the central flowline of Castaño Overa glacier, field measurements show an 8% decrease in surface velocity between March–May and May–September. Pléiades surface velocity shows the same behavior, with an 11% decrease in velocity. On Ventisquero Negro, ablation stakes (S.1 and S.2) show contrasting behavior, which, although not so contrasting, is also seen in the Pléiades surface velocity maps. S.2 shows a 5% decrease in velocity between March–May and May–September, while S.1 measurements indicate a 16% increase in velocity in the same period. For all four stakes, field and Pléiades measurements agree well within the error bars (Fig. 8a and b). The same agreement between field and satellite measurements is found for the ice flow directions (Table 5).
Fig. 8. (a) Comparison of Pléiades surface velocities from March to May 2012 with the displacement of ablation stakes measured in the field between March and May 2013. (b) Comparison of Pléiades surface velocity from May to June 2012 with the displacement of ablation stakes measured between May and September 2013.
Table 5. Comparison of velocity measurements at ablation stakes and derived from Pléiades imagery. The mean surface velocity difference is 0.003 m d−1 with a standard deviation of 0.03 m d−1. Mean direction difference is −1.1° with a standard deviation of 4.4°
Sub-seasonal changes in flow velocity
Although all periods have a similar glacier flow pattern, we found interesting sub-seasonal velocity changes on Monte Tronador glaciers. These changes are particularly evident in ablation areas where there is greater data coverage (Figs 6, 9 and 10). Gently sloping glaciers, which have a mean slope of <20° and an even longitudinal profile, show minor velocity fluctuations between March and June (unnamed glaciers 1, 2 and 3) or a minor decrease in mean surface velocity from April to June (Castaño Overa and Alerce glaciers) (Fig. 10a–e). Steeper glaciers, with mean slopes >20° and icefalls or changes in longitudinal profiles, show contrasting behaviors. Casa Pangue, Manso and Verde glaciers, which have debris-covered tongues below a steep icefall/ice cliff, show an acceleration from March/April to June (Fig. 10i–k). The acceleration is concentrated in their debris-covered tongues, and is stronger for Verde glacier. On the other hand, Frías, Norte and Blanco Chico glaciers, which do not descend below the rock cliff and do not have debris-covered tongues, show a decrease in surface velocity between March and June (Fig. 10f–h).
Fig. 9. Surface velocity profiles over Verde (a) and Frías (b) glaciers (see blue in Fig. 1).
Fig. 10. (a–k) Mean surface velocities for the central flowline of the ablation area of glaciers of Monte Tronador for the four different periods. Date format on horizontal axis is dd/mm (2012). The mean altitude of the ablation zone with ±1 standard deviation is given for each glacier. The light-color shaded area indicates the ±1 standard deviation of the velocity measurements. (l) Location of the profiles shown in (a–k). On the map, the colors of the profiles correspond to those shown in the surface velocity diagrams. The ablation areas of the glaciers are shown in dark gray and the accumulation areas are in light gray.
Discussion
Precision of the velocity measurements
Image resolution is one of the factors limiting the precision of the displacements retrieved from multitemporal analyses of satellite images (Reference Leprince, Barbot, Ayoub and AvouacLeprince and others, 2007). The use of images with sub-meter resolution, such as Pléiades, leads to a higher precision, implying that useful velocity fields can be derived from images acquired only a few days or weeks apart, even for slow-flowing glacier ice. This is a significant advantage because short-period changes of the glacier surface (e.g. changing in snow-cover distribution, ablation and crevasse morphology) are decreased, improving correlation and allowing the derivation of a dense and nearly comprehensive velocity field (Reference Riveros, Euillades, Euillades, Moreiras and BalbaraniRiveros and others, 2013). Monthly or weekly fluctuations in glacier velocities could now be observed from time series of optical satellite images, something only previously possible using radar images and for fast-flowing glaciers (Reference JoughinJoughin, 2002; Reference Scherler, Leprince and StreckerScherler and others, 2008; Reference Joughin, Smith, Shean and FloricioiuJoughin and others, 2014).
The main disadvantage of Pléiades images compared with other medium spatial resolution optical satellite images (e.g. ASTER or Landsat) is their cost (1 euro per km2 of image for European researchers subsidized by the ISIS program of the French Space Agency (CNES), http://www.isis-cnes.fr/). Nevertheless, being able to program the dates of image acquisition and the short revisit time of the two Pléiades satellites are particularly valuable advantages for monitoring special short-lived events and/or avoiding seasons with snow or high cloud cover.
The confidence intervals (Q(0.9)) of the velocity fields for March, April and May obtained here are between 0.03 and 0.06 m d−1 (11 and 22 m a−1 respectively). Errors are much larger for the June image pair (0.28 m d−1 equivalent to 102 m a−1) probably due to the increase in the snow-covered area. Based on the error analysis over bare ground not affected by cast shadows, we found that our Q(0.9) error estimate is conservative.
The precision of the surface displacements (roughly ±1 pixel or ±0.5 m at the 1σ confidence level, Q(0.683)) is not as high as reported in other similar studies where sub-pixel precisions are often reached (Reference Stumpf, Malet, Allemand and UlrichStumpf and others, 2014; Reference Lacroix, Berthier and MaquerhuaLacroix and others, 2015). Apart from the rough topography and the presence of cast shadows in the densely vegetated area surrounding Monte Tronador, a further reason for this relatively low precision is the large difference in incidence angles between the correlated images (the base-to-height ratio in Table 2 is always >0.34). Ideally, a pair of images with much smaller base-to-height ratio should be selected (or tasked) to reach sub-pixel precision (Reference DelacourtDelacourt and others, 2007). Such a pair of images was not available in the case of Monte Tronador glaciers.
Another possible source of error related to these relatively high base-to-height ratios is associated with glacier elevation differences between the date of the DEM and the dates of the images. This source of error does not occur on ice-free terrain and increases linearly with the elevation differences and base-to-height ratio (Reference BerthierBerthier and others, 2005; Reference Scherler, Leprince and StreckerScherler and others, 2008). There is no published information about recent elevation changes from glaciers of Monte Tronador. Instead, we compared our GNSS measurements (6000 data points) over Alerce and Ventisquero Negro glaciers on 16–18 May 2013 and the Pléiades DEM (21 April 2012). The maximum elevation changes were on Ventisquero Negro and reached −10 m in 390 days. Assuming that the surface lowering occurred at the same rate between the successive Pléiades images, we estimated a possible maximum error of 0.2–0.3 m, which is well within our error estimates based on the Q(0.9) over the ice-free terrain.
Reported errors of surface velocity measurements of other Patagonian glaciers (mainly outlet glaciers of the North and South Patagonian Icefields) obtained by cross-correlation or manual feature tracking of optical images are in the range 15–300 m a−1 (Reference Skvarca, Raup and De AngelisSkvarca and others, 2003; Reference Melkonian, Willis, Pritchard, Rivera, Bown and BernsteinMelkonian and others, 2013; Reference Sakakibara, Sugiyama, Sawagaki, Marinsek and SkvarcaSakakibara and others, 2013; Reference Sakakibara and SugiyamaSakakibara and Sugiyama, 2014). Even though the acquired conditions of Pléiades images were not optimal (high base-to-height ratio), the precision of our measurements (typically 10–20 m a−1) is higher or within the lower range of the velocity errors reported in previous studies.
Flow patterns
The surface velocity maps obtained are highly consistent over different periods and represent, with a spatial resolution of 16 m, the flow pattern of the glaciers of Monte Tronador. The highest velocities are located on the steeper parts of Verde and Manso glaciers and are likely associated with a high mass turnover on these sectors (Fig. 9a). In contrast, the lower reaches of the debris-covered tongues of Verde and Casa Pangue glaciers show the lowest velocities, <10 m a−1. This is in agreement with the conclusions of Reference BownBown (2004) and Reference Bown and RiveraBown and Rivera (2007), who considered that the lower tongue of Casa Pangue contains stagnant ice zones. The velocity maps obtained here also show many other interesting characteristics of these glaciers. For example, both Frías (Fig. 9b) and Ventisquero Negro (Fig. 8) glaciers exhibit a marked acceleration of ice flux at their fronts, likely related to the dry calving observed in the front of Frías and the lake calving at Ventisquero Negro. Although Castaño Overa and Alerce glaciers both terminate close to a cliff, they do not show acceleration close to the terminus. Thus, we suggest that they are unaffected by dry calving. This agrees with their lower retreat rate in recent years, compared with Frías or Manso/Ventisquero Negro glaciers.
Velocity maps show that the ice flow crosses ridges in the ice surface in some cases, suggesting that the use of hydrological routines in GIS to delimit ice divides from a DEM should be used with caution in areas with strong topographic differences (Fig. 7). In this regard, the flow patterns obtained from the surface velocity maps will provide a very valuable tool to aid in the delimitation of ice divides at Monte Tronador (inset in Fig. 7). Correct delimitation of the surface area of each glacier could help improve the derivation of glaciological parameters such as the accumulation–area ratio and the analysis of ice flux, mass turnover and mass balance of these ice masses.
Analysis of sub-seasonal changes in flow velocity
Changes in surface velocity from daily to seasonal periods have been successfully associated with change in subglacial water pressure (Reference Iken and BindschadlerIken and Bindschadler, 1986; Reference SugiyamaSugiyama and others, 2011). The water pressure under glaciers is controlled by characteristics of the glacier bed, the ice flow dynamics and the amount of water reaching the bed. Changes in surface melt rate, precipitation and/or water depth for calving glaciers cause complex spatial and temporal flow variations (Reference Sugiyama and GudmundssonSugiyama and Gudmundsson, 2004). Empirical results show that glaciers accelerate as water pressure increases (Reference Iken and BindschadlerIken and Bindschadler, 1986; Reference Sugiyama and GudmundssonSugiyama and Gudmundsson, 2004), with an intra-annual cycle, controlled by ablation, meltwater generation and the efficiency of the subglacial hydrological system, that reaches a maximum at the end of the ablation season (Reference Nienow, Sharp and WillisNienow and others, 1998). At smaller timescales, particular events such as heavy rainstorms or the sudden drainage of supraglacial lakes can cause surface speed to increase for a few hours (Reference Naruse, Fukami and AniyaNaruse and others, 1992; Reference Sugiyama and GudmundssonSugiyama and Gudmundsson, 2004).
Unfortunately, at present there is no information about the water runoff or subglacial water pressure variations of glaciers of Monte Tronador. Thus, we developed a simple model to identify periods when more meltwater or rain could reach the subglacial system (see Appendix). Using a positive degree-day approach and data from a nearby meteorological station, the model calculates the meltwater generated at the glacier surface and discriminates between solid and liquid precipitation, allowing the analysis of the two main sources of water (meltwater and rain) to the subglacial system.
After calibrating the model with mass-balance data from stakes 1 and 4 for 2013 (see Appendix), we ran the model for 1700 and 1050 m a.s.l. during the time covered by the Pléiades scenes (Fig. 11). The highest elevation approximately corresponds to the ablation areas of Frías, Norte and Blanco Chico glaciers, whereas the lower elevation represents the altitude of the debris-covered tongues of Manso, Casa Pangue and Verde glaciers (Fig. 11a and b).
Fig. 11. Results of the mass-balance model for two different elevations at Monte Tronador, covering the same time span as the Pléiades pairs analyzed in this study. (a) Mass-balance model results at 1700 m a.s.l., the mean elevation of the ablation area of highly sloping, upper clean ice glaciers. (b) Mass-balance model results at 1050 m a.s.l., the mean elevation of the ablation area of low-slope, debris-covered, lower-elevation glaciers.
The model shows that, in the ablation area of the higher-elevation clean-ice glaciers (1700 m a.s.l.), the amount of liquid water decreased from March to June, following the expected seasonal drop in temperature. According to the model, at this elevation, in June there was almost no meltwater generation or rain because the low temperatures inhibited melting and promoted solid precipitation (Fig. 11a). On the other hand, in the ablation area of the lower-elevation debris-covered glaciers (1050 m a.s.l.), the decrease in meltwater production from March to June due to decreasing temperatures was balanced by the occurrence of rainstorms in May and, especially, in June (Fig. 11b). The increase of available water from these rainstorms could explain the increase in surface velocity observed in May and June on these lower glacier tongues. In comparison to the other debris-covered glaciers, Verde glacier shows a greater increase in surface velocity, which is related to acceleration in the central area of its tongue. This area contains a clearly visible rock/debris avalanche deposit that fell from the southwest margin of the bedrock cliff and has been subsequently transported on the surface of the glacier. It is possible that ice thickness is greater here (indeed, the longitudinal elevation profile of Figure 9a shows a rise in this area) due to the thermal insulation effect of the avalanche deposits. Since subglacial motion depends on the basal shear stress, which depends on the ice thickness and surface slope (Reference Cuffey and PatersonCuffey and Paterson, 2010), a higher ice thickness appears to be a plausible explanation for the higher surface velocity in this area.
Despite these findings, it is clear that further research is needed to improve understanding of the relationship between changes in the hydrological system of these glaciers and their impact on basal motion. For example, several studies (e.g. Reference Nienow, Sharp and WillisNienow and others, 1998; Reference Sugiyama and GudmundssonSugiyama and Gudmundsson, 2004) have shown that the development of the hydrological system of a glacier over the ablation season depends on the transport time for meltwater to reach the bed. Those authors also indicate that the type of channel connection at the glacier bed can have major implications for the seasonal to hourly glacier surface velocity variations. In our study, velocity measurements only started near the end of the ablation season (March), and ideally the velocity fields should be available over the entire summer to investigate differences in dynamic responses of the glaciers to surface melting. Continuous GNSS measurements of the surface velocity at selected locations on the glaciers would provide further insight into the relationships between short-term velocity variations and liquid water availability.
Conclusions
In this study, we generated accurate, high-resolution seasonal surface velocity maps of several glaciers by applying cross-correlation techniques to sequential Pléiades imagery. We measured ice surface displacements on glaciers at Monte Tronador over periods as short as 19 days, with a precision of ∼0.58 m (i.e. 11 m a−1). These are the first surface velocity measurements and maps for glaciers of Monte Tronador, and the first such maps derived from Pléiades imagery. Our results show that these glaciers flow radially from the summit, with maximum surface speeds of ∼390 m a−1 associated with steep icefalls. The lower reaches of the debris-covered tongues of Verde and Casa Pangue glaciers are very low and could be considered stagnant ice, whereas another debris-covered tongue (Ventisquero Negro) shows acceleration at its lake-terminating calving front.
The glaciers showed different patterns of behavior over the period March to June 2012, which were tentatively related to differences in the amount of meltwater and/or rain reaching the glacier bed. Low-elevation debris-covered glacier tongues (Manso, Verde and Casa Pangue) showed an increase of surface velocity at the beginning of the accumulation season (May–June), probably related to the percolating water from rainfall events. At higher elevations, the clean-ice tongues of Frías, Blanco and Norte glaciers decreased their velocities over the same time span, probably in response to a decrease in temperatures and meltwater production. Gently sloping glaciers (Castaño Overa, Alerce and three unnamed glaciers in Chile) show little surface velocity fluctuation during our study period.
Such detailed glacier surface velocity maps have several glaciological applications and constitute helpful tools to demarcate the boundaries between adjacent glaciers in complex terrain. However, a more continuous monitoring of glacier surface velocity and water pressure measurements is needed to carry out detailed analyses of the sub-seasonal changes in surface velocities of these glaciers.
Acknowledgements
We acknowledge the support from Secretaría de Ambiente y Desarrollo Sustentable de Argentina, Agencia de Promoción Científica (projects PICT 2010-1438) and Consejo Nacional de Investigaciones Científicas y Técnicas, Argentina (CONICET), for funding. Etienne Berthier acknowledges support from the French Space Agency (CNES) through the TOSCA program. Pléiades images were provided at no cost by Airbus Defense and Space through the Pléiades User Group. Administración de Parques Nacionales kindly provided permission and logistical assistance to work at Monte Tronador inside Parque Nacional Nahuel Huapi. The mass-balance and stake position measurements used in this study could not have been conducted without the valuable field knowledge and collaboration provided by Mauricio Cadillo, Ernesto Corvalán, Jorge Luis Lambert and Juan Pablo Scarpa. The fieldwork would have been much harder without the hospitality of Nicolas Betinelli and the staff of Refugio Otto Meiling. The constructive comments of two anonymous reviewers, the Scientific Editor (David Rippin) and Brian Luckman helped to improve the manuscript. We acknowledge the assistance of the International Glaciological Society (IGS) in the publication of our paper.
Appendix
Here we describe the temperature index mass-balance model used to identify periods when meltwater or rain could reach the subglacial system of Monte Tronador glaciers. Input data are from Bariloche Airport weather station (GHCND:AR000087765; 41.05° S, 71.07° W; 840 m a.s.l.), located ∼60 km east of Monte Tronador.
The model is based on the positive degree-day method (Reference HockHock, 1999) in combination with a simple parameterization to discriminate between rain and snow (Reference OerlemansOerlemans, 2001). Melt rate M (m w.e. d−1) at each selected elevation is obtained by the basic formulation of a temperature index approach:
where DDF is the degree-day factor (m d−1 °C−1), different for snow, ice and debris-covered ice surface, z is the selected elevation (m) and T (z) is the daily mean air temperature at each specific elevation:
where z 0 = 840 m a.s.l., the altitude of Bariloche station, l is the temperature lapse rate and T 0 is the mean daily air temperature at Bariloche station (°C).
Precipitation P (z) at each selected elevation is obtained by assuming a linear increase with elevation z:
where s is the lapse rate of precipitation and Pz 0 is the total daily precipitation at Bariloche station (m). We use a fixed air temperature threshold (T t) to determine whether precipitation falls as snow or rain.
To calibrate the model we used data from ablation stakes S.1 and S.4 installed on Monte Tronador in March 2013. S.1 was measured three times and S.4 was measured seven times between 8 March and 23 September 2013 (Fig. 12). Stake S.1 is located in the debris-covered tongue of Manso glacier at 1050 m a.s.l. Stake S.4 is located close to the ice divide between Castaño Overa and Alerce glaciers at 2020 m a.s.l. (see Fig. 1). We used mean daily air temperatures and total daily precipitation from Bariloche and calibrated the model interactively by gradually adjusting the different parameters (DDFsnow/ice/debris, l, s or T t). The best combination of parameters is presented in Table 6. The RMSE between the mass-balance measurements and the model results for stakes S.1 and S.4 were 0.25 and 0.15 m w.e., respectively. The differences between the model results and our measured data could be because this simple approach does not take into consideration the spatial variation associated with slope, aspect or wind direction in ablation and/or accumulation processes. Nevertheless, despite the simplistic nature of this approach, the good fit of the modeled values suggests that these equations can capture the main mass-balance variations observed on the glaciers of Monte Tronador (Fig. 12).
Fig. 12. Measured and modeled accumulated mass balances at stakes S.4 (a) and S.1 (b).
Table 6. Best combination of parameters to reproduce the accumulated mass balance measured between March and September 2013 at stakes S.1 and S.4
