Introduction
Mountain glaciers are highly sensitive to variations in temperature and precipitation (Reference OerlemansOerlemans, 2005) and thus represent a good indicator of climate changes over recent centuries (e.g. Reference Oerlemans and FortuinOerlemans and Fortuin, 1992; Reference Haeberli, Maisch and PaulHaeberli and others, 2002). Furthermore, volume changes in glaciers and small ice caps, excluding the Antarctic and Greenland ice sheets, are projected to make the largest contribution, after thermal expansion of the world’s oceans, to sea-level rise in the coming century (Reference MeierMeier and others, 2007). Many of the world’s glaciers are in remote areas, meaning that land-based methods of measuring ice volume and velocity changes are expensive and labour-intensive. Remote-sensing technologies have offered a solution to this problem through the use of satellite radar interferometry and feature tracking using both optical and radar imagery (e.g. Reference KääbKääb, 2002, Reference Kääb2005; Reference BerthierBerthier and others, 2005; Reference Rignot and KanagaratnamRignot and Kanagaratnam, 2006; Reference Luckman, Quincey and BevanLuckman and others, 2007; Reference Quincey and GlasserQuincey and Glasser, 2009). The efficacy of using feature tracking over short time windows to obtain a complete and detailed map of ice velocity on mountain glaciers has been demonstrated in recent studies (e.g. Reference KääbKääb, 2002, Reference Kääb2005; Reference BerthierBerthier and others, 2005; Reference Scherler, Leprince and StreckerScherler and others, 2008). Precise orthorectification of optical imagery appears to be the key to obtaining whole-glacier coverage of ice velocities. Reference Leprince, Barbot, Ayoub and AvouacLeprince and others (2007a,Reference Leprince, Ayoub, Klinger and Avouacb, Reference Leprince, Berthier, Ayoub, Delacourt and Avouac2008) and Reference Scherler, Leprince and StreckerScherler and others (2008) have obtained sub-pixel orthorectification and correlation accuracy and demonstrated its applicability in measuring glacier velocities.
Large changes in mountain glacier dynamics are usually only apparent over periods of decades. As a result there are relatively few studies which examine interannual variations in ice velocity (Reference Lliboutry and ReynaudLliboutry and Reynaud, 1981; Reference Span, Kuhn and SchneiderSpan and others, 1997; Reference Vincent, Vallon, Reynaud and MeurVincent and others, 2000, Reference Vincent, Soruco, Six and Le Meur2009; Reference Span and KuhnSpan and Kuhn, 2003). These studies provide details of the decadalscale relationship between mass balance, ice velocity and, at times, glacier advance and retreat. While the time span of these studies is long, the spatial resolution of measurement is necessarily low.
We present here the application of a high-precision orthorectification, co-registration and correlation technique (Reference Leprince, Ayoub, Klinger and AvouacLeprince and others, 2007b) to obtain horizontal velocities over temperate glaciers. We apply the technique to fast-flowing glaciers of the central Southern Alps of New Zealand and compare these observations to GPS ground-based velocity measurements. The maritime climate and high mass turnover of these glaciers means that their velocities are very high for alpine glaciers (i.e. up to ∼5 m d−1). Hence it is possible to obtain velocities from all but the slowest-moving parts of the accumulation and ablation areas by comparing images separated by relatively short periods. Furthermore, these glaciers are amongst the few worldwide that have been through a retreat/advance phase during the current period of Earth warming, and this change has occurred during the relatively short period of satellite observation. Changes in glacier velocities through time can provide insights into the differing response of mountain glaciers to climatic forcing by illustrating how glacier basal velocities vary during a retreat/advance cycle.
Study Area
The Southern Alps of New Zealand present a ∼600 km long and ∼2.5 km high barrier to the prevailing westerly winds, resulting in high precipitation (e.g. Reference Henderson and ThompsonHenderson and Thompson, 1999), cool temperatures and temperate glacierization. Approximately 1158 km2 of glacier surface area and ∼60 km3 of ice volume is concentrated in the central part of the Southern Alps (Reference ChinnChinn, 2001). A few large glaciers contain most of the Southern Alps ice volume (Reference ChinnChinn, 2001), for example Tasman, Murchison, Fox and Franz Josef Glaciers, with 29%, 7%, 6% and 6%, respectively, of the total volume. This study focuses on Franz Josef, Fox and Tasman Glaciers and the area around them (Fig. 1).
Since they were first measured close to their final Little Ice Age maxima in the late 19th century, glaciers in the Southern Alps have lost significant volumes of ice (Reference ChinnChinn, 1996; Reference Hoelzle, Chinn, Stumm, Paul and HaeberliHoelzle and others, 2007; Reference SchaeferSchaefer and others, 2009). This continues today. A large fraction of this ice loss is caused by the formation and growth of proglacial lakes, the final result of decades of downwasting (Reference Kirkbride and WarrenKirkbride and Warren, 1999). The most comprehensive record of glacier length fluctuation has been established for Franz Josef Glacier, although good records also exist for Fox and Tasman Glaciers (Fig. 2).
Within this loss of ice mass, there has been a divergence of individual glacier response, with adjacent glaciers exhibiting opposite terminus responses. While the glaciers have all experienced a major retreat since the late 19th century, the timing of retreat has varied by decades. Indeed, the largest glacier in the Southern Alps, Tasman Glacier, has steadily thinned but did not retreat until the start of proglacial lake growth in 1983, while more responsive glaciers, such as Fox and Franz Josef Glaciers, have experienced more complex advance/retreat histories. Figure 2a shows that Fox and Franz Josef Glaciers have retreated substantially since the late 1930s, but includes a substantial advance since the mid-1980s. In the last decade, these glaciers have experienced a period of retreat followed by an advance (300–400 m variations; Fig. 2b). Their advance in recent years is an unusual phenomenon that provides a unique opportunity to examine the dynamics of glacier advance and retreat.
Methods
Remotely sensed displacement measurements
We quantify glacier surface displacement by performing orthorectification, co-registration and correlation of optical imagery, using the software, COSI-Corr (Co-registration of Optically Sensed Images and Correlation), described in detail by Reference Leprince, Barbot, Ayoub and AvouacLeprince and others (2007a) and Reference Scherler, Leprince and StreckerScherler and others (2008). Following the standard COSI-Corr procedure, a digital elevation model (DEM) was used as the global ground truth for the study of ASTER (Advanced Spaceborne Thermal Emission and Reflection Radiometer) imagery because no global georeference was available. We used a 50 m spatial resolution DEM provided by Terralink. The January 2006 acquisition is used as the first master image, and tie points are selected on stable features between the raw ASTER image and a shaded version of the DEM. These tie points are optimized using sub-pixel correlation and converted into ground-control points (GCPs) for precise orthorectification. Over 30 GCPs were used and the average mis-registration between this ASTER ortho-image and the DEM is 35 cm, with a standard deviation of 8.6 m. This very accurate co-registration between the DEM and the imagery ensures accurate modelling of topography effects in the ortho-image. By construction, the absolute georeferencing accuracy of the image is that of the DEM. Tie points are then selected on stable features between this first master orthoimage and subsequent raw ASTER images. Tie points for each image are optimized using sub-pixel correlation and converted to GCP for orthorectification. Eight GCPs are generated between each pair of images. On average, mis-registration is <2 cm, with a standard deviation of <3.5 m. This large deviation is produced by unrecorded attitude jitter of the satellite platform.
For this study we only selected nadir-looking ASTER acquisitions (Table 1). The velocity estimates presented are therefore free from any potential biases, which could be due to elevation model errors or from potential ice-thickness changes between acquisitions. The VNIR (visible/near-infrared) bands 3N at 15 m ground resolution were used. Horizontal displacements were measured from sub-pixel correlation using a multi-scale approach where the smallest correlation window size is 32 × 32 pixels, sliding every pixel. A destriping of the displacement field was performed to remove artifacts generated by undulations of the observing platforms (Reference Leprince, Ayoub, Klinger and AvouacLeprince and others, 2007b; Reference Scherler, Leprince and StreckerScherler and others, 2008; Reference Teshima and IwasakiTeshima and Iwasaki, 2008). Finally, reported uncertainties represent the standard deviation of the measurements in each north–south and east–west component, over areas of stable ground, as explained in detail by Reference Scherler, Leprince and StreckerScherler and others (2008).
Ground-based velocity measurements
Ice velocity was measured at Franz Josef Glacier by repeat survey of a network of GPS stakes between November 2001 and April 2002 (location of GPS stations shown in Fig. 1 and dates of measurement given in Table 2). The stakes are arranged in one longitudinal line and two transverse lines, although there is a large gap in longitudinal stake coverage due to extensive crevassing. Most velocity measurements were acquired after image acquisition (i.e. from February to April 2002; Table 2); however, at four stakes the velocity measurements started before image acquisition and overlapped with ASTER-derived displacement measurements.
The surveying was undertaken using Trimble 4700 GPS receivers. The GPS data were processed using Trimble software, which typically reported horizontal errors in position in the range 0.01–0.02 m. The maximum horizontal error in stake displacement between two measurements (0.2 m) includes a GPS error of 0.04 m; the remainder is antenna positioning error.
Results
In Figures 3 and 4 we present the Euclidean norm of the displacements measured in the north–south and east–west directions over a 16 day window in summer 2002 and 2006, respectively. For both measurements, large displacements are observed, indicating large velocities integrated over a 16 day time window, especially for Franz Josef and Fox Glaciers (i.e. up to ∼5 m d−1). High marginal velocities suggest that basal sliding dominates ice motion in the ablation zones of the steeper glaciers. It is also interesting to note that the images highlight very similar displacement patterns between the two periods.
We now compare the ASTER-inferred velocities with the ground-based velocity measurements on Franz Josef Glacier (Table 2; Fig. 5). The velocities derived from the correlation measurements were linearly interpolated at the location of the GPS stations. Given the high resolution of the ASTER images (15 m), the uncertainties related to the interpolation should be negligible. We observe a systematic overestimation of the velocities in the satellite data (Fig. 5a). While this mismatch is expected since many of the stake measurements were made later in the ablation season, when sliding velocities are expected to be lower, the overestimation appears to be large where velocities are lower (i.e. <0.5 m d−1), regardless of being in the accumulation or ablation area or of the date of acquisition. As noted by Reference Leprince, Ayoub, Klinger and AvouacLeprince and others (2007b) and Reference Necsoiu, Leprince, Hooper, Dinwiddie, McGinnis and WalterNecsoiu and others (2009), the noise on each east–west (EW) and north–south (NS) component of the correlation can be modelled as additive white Gaussian noise. Given that x and y are two independent random variables following a normal distribution with variance σ 2, the magnitude given by follows a Rayleigh distribution with mean (Reference MeikleMeikle, 2008). Computing the Euclidean norm can, in turn, lead to an overestimate close to μ. We obtain σ = 2.6 m, hence μ = 3.2 m. Given the uncertainty obtained here, the mismatch between ground-based and satellite-derived velocities could be related to this for velocities less than 0.3 m d−1. Furthermore, the discrepancies can be very large when the displacements are very small (i.e. velocities less than 0.1 m d−1 over the 16 day period). To preserve the actual noise properties (i.e. zero mean and additive white Gaussian noise) and avoid potential overestimation of the velocities, it is therefore more appropriate to compare velocities in the north–south and east–west directions. The results are shown in Figure 5b and c. They are similar to those for the computed norm, meaning that they show discrepancies for small velocities. A common bias generally observed in glacier velocities derived from satellite cross-correlation is due to glacier thickness changes between image acquisitions. Parallax errors then occur in the epipolar plane defined by the viewing geometry of each image (Reference BerthierBerthier and others, 2005; Reference Scherler, Leprince and StreckerScherler and others, 2008). For our study, we computed that the epipolar directions, projected on the ground, had an orientation of N79.1° W, which is mostly east–west oriented. Hence, if glacier thickness changes had biased our estimations in 2002 (with respect to the GPS measurements), the EW component of the velocity field would have suffered most of this bias. However, as observed in Figure 5, the EW component shows little bias. It shows a zero-centred dispersion of the measurements below velocities of 0.2 m d−1, which corresponds to average displacements of 3.2 m over the observation period and is consistent with the noise standard deviation of 2.2 m estimated from Figure 3. Indeed, because of the very close nadir acquisitions (Table 1), the horizontal velocities derived are insensitive to glacier thickness changes. For instance, the glacier thickness would need to have changed by >12 km during the 16 day observation period to induce a bias of >1 m, or 6 cm d−1 (using equation (1) of Reference Scherler, Leprince and StreckerScherler and others, 2008). The bias observed in the NS component is therefore not due to glacier thickness changes. However, close inspection reveals a slight attitude residual from the spacecraft, which may not have been completely corrected during our a posteriori attitude-removal procedure. We therefore explain the bias in the NS component by small pitch residuals of the spacecraft for which we were unable to correct.
The ice velocities measured on Tasman Glacier are similar to those obtained using some terrestrial and photogram-metric measurements (Reference Kirkbride and WarrenKirkbride and Warren, 1999) and ASTER-derived measurements (Reference KääbKääb, 2002; Reference Quincey and GlasserQuincey and Glasser, 2009), although no GPS-based measurements were acquired there. For the lowermost 2–3 km of the glacier terminus, the velocities observed from ASTER data seem to be similar, within the measurement accuracy, to those observed for 1971–86 and 1957–71. Terrestrial velocity measurements from May to December 1986 at a region ∼8 km up-glacier of the front gave 0.24 m d−1, which agrees well with both 2002 and 2006 measurements.
We continue our analysis of the measurements by comparing velocities along flowlines for the three glaciers (Franz Josef, Fox and Tasman Glaciers). The flowlines are shown in Figure 4 and the extracted velocities in Figure 6. There is no detectable change in velocities for Tasman Glacier, which has been in a state of continual retreat (Fig. 6a). While small changes in velocities are observed along most of Franz Josef and Fox Glaciers, large changes in velocities (up to a threefold increase in both instances) are observed at the snouts of the glaciers during the advancing phase (Fig. 6b and c). This magnitude of velocity change must, however, be treated with caution because the variation is observed for velocities around and lower than 0.5 m d−1. Finally, it is also interesting to note that the maximum velocities for Franz Josef and Fox Glaciers occur below the location of the mean estimated equilibrium-line altitude (ELA) at ∼1860 m a.s.l. (Reference Anderson, Lawson, Owens and GoodsellAnderson and others, 2006).
Discussion
A significant change in integrated velocity between 2002 and 2006 has been identified in the lower ablation zone of Franz Josef and Fox Glaciers by repeated calculation of glacier integrated velocities using ASTER imagery. One possible explanation of this is that the observed changes in velocity are a result of taking short snapshots of velocity within a velocity field that is highly variable in time. Previous work (Reference AndersonAnderson, 2003) has identified velocity variations at sub-monthly and, to a lesser extent, seasonal scales. As both sets of measurements in the present study were taken at the same time of year (midsummer) we can neglect seasonal velocity variability as a reason for the change in velocity observed.
Many mountain glaciers exhibit short-term (sub-monthly) velocity variations due to variations in sliding velocities (e.g. Reference Iken and BindschadlerIken and Bindschadler, 1986). Such variations have been observed at Franz Josef Glacier throughout the year in response to large inputs of water, from melt or rainfall. In both 2002 and 2006 there was significant rainfall in early January, and the totals for the 6 week period from the start of January to mid-February (covering the period of both measurements and 3–4 weeks beforehand) were similar (Fig. 7). The melt rates on the lower part of the glacier, measured as part of the programme described by Reference Anderson, Lawson, Owens and GoodsellAnderson and others (2006), were 20% lower in 2006 than 2002. Glacier velocity is not simply related to the quantity of melt and rainwater but also depends in a complex way on the effective pressure of water at the bed. Generally, high water pressures occur when large amounts of water enter a poorly developed or inefficient subglacial drainage system. The large amounts of water available in January–February 2002 and 2006 do not suggest this was the case in either of the measurement periods.
The sliding velocity we observe cannot be reconciled with an empirical law where sliding velocity on a hard bed is proportional to basal stresses and inversely proportional to the effective pressure (Reference PatersonPaterson, 1994 and references therein):
where u b is the basal sliding, C s is a sliding parameter, τ b is the basal shear stress, N is the effective pressure (equal to p i − p w, where p i and p w are the ice overburden and water pressure, respectively), and p and q are constants. This law is commonly used in models of temperate glaciers (e.g. Reference BindschadlerBindschadler, 1983; Van der Veen, 1987; Reference PattynPattyn, 1996; Reference Le Meur and VincentLe Meur and Vincent, 2003) and implies a basal drag that increases without bound with increasing sliding velocity or effective pressure. Assuming a vertically integrated gravity-driven shear stress, τ b can be approximated as being proportional to the ice thickness and the slope at the surface of the ice. The increased thickness alone (average 67 m thickness along measured profile in 2002, average 121 m thickness in 2006) would increase the velocity by a large factor, due to internal deformation. The observed changes are much less than this. Alternatively, a physically based sliding law, in the presence of cavitation, initially proposed by Reference LliboutryLliboutry (1968, Reference Lliboutry1987) and Reference FowlerFowler (1986) and recently revived by Reference SchoofSchoof (2005) and Reference Gagliardini, Cohen, Råback and ZwingerGagliardini and others (2007), expresses the relationship between basal friction and sliding velocity as:
This relationship has been described as a function of three independent parameters that only depend on the bedrock geometry (Reference Gagliardini, Cohen, Råback and ZwingerGagliardini and others, 2007). Unlike the empirical law described above, this phenomenological law imposes an upper limit to basal shear stresses and a decrease at low effective pressure and high velocity, as discussed in detail by Reference Gagliardini, Cohen, Råback and ZwingerGagliardini and others (2007). Given the very high velocity and likely low effective pressure observed here, basal shear stress may be relatively insensitive to changes in ice thickness. It would, however, be more sensitive to changes in local water discharge. The local water budget measurements we have acquired here unfortunately do not allow us to distinguish between these behaviours. Piezometers should be used to monitor changes in effective pressure.
Another, more likely, possibility is that the observed changes in velocity are caused by longer-term changes in glacier geometry driven by mass-balance changes. In this case both the dynamic response and climatic forcing of these glaciers must be considered in order to understand the contrasting behaviour between Fox/Franz Josef and Tasman Glaciers. First, the differences in dynamics between these glaciers can be partially captured by assessing the response time. Franz Josef Glacier and the neighbouring Fox Glacier have exceedingly short response times estimated at ∼7 years (estimated using a characteristic ice thickness of ∼150 m and an ablation rate at the terminus of ∼20 m a−1; Reference Jóhannesson, Raymond and WaddingtonJóhannesson and others, 1989; Reference Anderson, Lawson, Owens and GoodsellAnderson and others, 2006). Tasman Glacier’s response time is hard to estimate using this method because the debris covering a large part of the ablation area reduces the terminus ablation rate by ∼90% (Reference KirkbrideKirkbride, 1989; Reference Purdie and FitzharrisPurdie and Fitzharris, 1999). Response time estimates between 20 and 200 years can be obtained using an ice thickness of ∼500 m and ablation values between 2.5 m a−1 under thick debris cover (Reference Kirkbride and WarrenKirkbride and Warren, 1999) and ∼22 m a−1 on clean ice near the terminus (Reference RöhlRöhl, 2008). Tasman Glacier’s response is also complicated by the growth of a proglacial lake (Reference Purdie and FitzharrisPurdie and Fitzharris, 1999) which, to some extent, has decoupled ice loss at the terminus from climatic influence. Second, the climatic forcing for these three glaciers is essentially the same (Purdie and others, in press), although that forcing is expressed differently for each glacier through its annual mass balance. While there have been many point measurements of mass balance (Reference ChinnChinn, 1969; Reference Anderson, Lawson, Owens and GoodsellAnderson and others, 2006; Reference Purdie, Brook and FullerPurdie and others, 2008), there are no complete mass-balance records for this area. A mass-balance record reconstructed from lowland climatic data (Reference Anderson, Lawson, Owens and GoodsellAnderson and others, 2006) for Franz Josef Glacier (Fig. 8) shows that in recent decades there have been three periods of positive balance in each of the 1980s, 1990s and 2000s, which can be matched with the length record shown in Figure 2. Since 1984, Franz Josef Glacier has been in a general phase of advance, although this advance has been punctuated by several short standstill and retreat events (Fig. 2b). One of these retreat phases lasted from 1999 to 2005, during which the first velocity snapshot was taken, while from 2005 to mid-2008 the glacier was advancing, when the second velocity snapshot was taken. Coincidentally, the glacier was almost exactly the same length in 2002 and 2006. However, longitudinal profiles of glacier elevation collected at the snout using GPS in 2002 and 2006 show an increase in thickness of ∼80% between 2002 and 2006 (Fig. 9). This thickening is associated with an increase in surface slope near the terminus, in turn leading to increased velocities. Our interpretation is, therefore, that the velocity change is the result of the glacier tongue thickening and slope increase as it changes from a retreating to an advancing phase.
Conclusions
Precise orthorectification and sub-pixel correlation of ASTER satellite imagery, used to derive the horizontal velocity field for the large glaciers in the Southern Alps of New Zealand, highlight very high ice velocity on large parts of these glaciers in 2002 and 2006. These high velocities allow short time windows between repeat imagery to be used to gather glacier velocities from large parts of the surface. The image processing enables velocities to be captured even in accumulation areas, where velocities are lowest and surface features ephemeral. However, a detailed comparison with GPS ground measurements shows that discrepancies can appear when the displacements reach the uncertainty of the satellite-derived displacement measurements (i.e. 2–4 m). Here, we attribute these differences to a slight attitude residual from the spacecraft, for which we could not correct.
A comparison between velocities in 2002 and 2006 highlights velocity variations that are most pronounced on the extremely sensitive and responsive ablation parts of the steep Fox and Franz Josef Glaciers. These velocity changes are related to glacier response during a retreat/advance cycle between 2002 and 2006. The increases in velocity observed at the termini of these responsive glaciers are the result of positive mass balance, leading to thickening, increase in slope and advance. No change in velocity is observed for larger slow-response, low-angled and debris-covered glaciers, such as Tasman Glacier. It is difficult to reconcile our observations with existing theories on basal sliding.
Ice-velocity changes occur on an interannual scale in response to mass-balance forcing. Mountain glaciers grow thicker and steeper during periods of positive mass balance, resulting in increased ice flux, and, if the increase is large enough, this causes glacier advance. There are relatively few studies that document these changes, as a result of sparsity of measurements and general global recession of glaciers in the 20th century. In 100 years of records, Reference Span, Kuhn and SchneiderSpan and others (1997) report three phases of velocity increase on Hinter-eisferner, Austria, but only one (in 1992) led to an advance. Our observations for Franz Josef and Fox Glaciers agree with those of Reference Vincent, Soruco, Six and Le MeurVincent and others (2009) who analysed 32 years of velocity data at Glacier d’Argentière, France, and observed the coincidence of increasing velocity and glacier advance. However, they differ from the findings of Reference Vincent, Vallon, Reynaud and MeurVincent and others (2000), of an increase in velocity at Glacier de Saint-Sorlin, France, during glacial retreat and ice thinning.
Acknowledgements
The surface elevation data on Franz Josef Glacier were collected by I. Owens. Terminus position data from Fox Glacier were collected by H. Purdie, and the Tasman Lake growth data were recorded by T. Chinn. S.L. was partially supported by US National Science Foundation (NSF) grant EAR-0636097 and the Gordon and Betty Moore Foundation. We thank B. Brock, M. Gudmundsson and an anonymous reviewer for their constructive reviews.