[go: up one dir, main page]

Academia.eduAcademia.edu
Remote Sensing of Environment 112 (2008) 3806–3819 Contents lists available at ScienceDirect Remote Sensing of Environment j o u r n a l h o m e p a g e : w w w. e l s e v i e r. c o m / l o c a t e / r s e Glacier-surface velocities in alpine terrain from optical satellite imagery—Accuracy improvement and quality assessment Dirk Scherler a,⁎, Sébastien Leprince b, Manfred R. Strecker c a b c Institut für Geowissenschaften, Universität Potsdam, 14415 Potsdam, Germany Electrical Engineering Department, California Institute of Technology, Pasadena, USA Institut für Geowissenschaften, Universität Potsdam, 14415 Potsdam, Germany A R T I C L E I N F O Article history: Received 11 March 2008 Received in revised form 19 May 2008 Accepted 31 May 2008 Keywords: Mountain glaciers Glacier velocity Himalaya Optical imagery Orthorectification Co-registration COSI-Corr ASTER SRTM-error A B S T R A C T The worldwide retreat of mountain glaciers has important consequences for the water, food, and power supply of large and densely populated areas in South and Central Asia. Successful mitigation of the hydrological impacts on societies as well as assessing glacier-related hazards require large-scale monitoring of glacier dynamics. However, detailed glaciological data from the Asian highlands are lacking, due to its size and difficult accessibility. We have applied a novel technique for precise orthorectification, co-registration, and sub-pixel correlation of Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) satellite imagery to derive surface velocities of Himalayan glaciers. Our approach allows for the correction of offsets due to attitude effects and sensor distortions, as well as elevation errors if a digital elevation model (DEM) from the Shuttle Radar Topography Mission (SRTM) was used for orthorectification. After postprocessing, the error on the displacements is on the order of 2–4 m per correlation. Translated into annual velocities, this error is reduced (increased) when the correlated images are more (less) than a year apart. Through application of a filtering procedure and several quality tests, the consistency of the results is validated to provide confidence in the remotely sensed velocity measurements, despite the lack of ground control. This novel approach allows fast, easy, and economically viable acquisition of detailed glaciological data in areas of difficult access and provides a means for large-scale monitoring of glaciers in high mountainous terrain. © 2008 Elsevier Inc. All rights reserved. 1. Introduction The global warming of climate has continued to cause the retreat of glaciers in many mountainous regions, and even the most optimistic scenarios for future temperature change involve pronounced glacier retreat over many decades to come (e.g., Oerlemans, 1994; IPCC, 2007a). This has important consequences for the global hydrological cycle, particularly in climatic threshold areas characterized by water stress. For example, the water, food, and power supply of densely populated regions in South and Central Asia are to a large degree dependent on snow and glacier melt water (Karim & Veizer, 2002; Winiger et al., 2005; IPCC, 2007b). Successful mitigation of the climate-related hydrological changes and their impacts on society therefore poses a pressing challenge, which calls for large-scale monitoring of glaciers and a better understanding of their dynamics (e.g., Haeberli et al., 2000, 2007; Kargel et al., 2005). Due to the large extent and difficult accessibility of high mountainous terrain, especially in Asian orogens, remote-sensing techniques provide an ⁎ Corresponding author. E-mail address: dirk@geo.uni-potsdam.de (D. Scherler). 0034-4257/$ – see front matter © 2008 Elsevier Inc. All rights reserved. doi:10.1016/j.rse.2008.05.018 efficient way to collect data in disparate regions. For example, satellite images have been used to track changes in glacier geometry (e.g., Paul et al., 2002; Khalsa et al., 2004; Aizen et al., 2007); analyze and monitor supraglacial lakes (Wessels et al., 2002); determine the equilibrium line altitude (Rabatel et al., 2005), and estimate annual mass balances of glaciers (Berthier et al., 2007). Remote-sensing tools can also be efficiently used to determine the ice velocity of a glacier, which is a particularly crucial variable because it determines ice discharge (e.g., Scambos et al., 1992; Goldstein et al., 1993; Joughin et al., 2004, Rignot & Kanagaratnam, 2006). Although glacier-surface velocities can be measured directly on the glacier with high accuracy at arbitrary spatial and temporal resolutions (e.g., Hubbard & Glasser, 2005), observations over long periods involve frequent revisits of the survey points, which can only be located on the accessible parts of a glacier. Therefore, field measurements commonly result in very sparse spatial coverage. In contrast, remote sensing-based measurements provide the opportunity to achieve large and possibly complete spatial coverage, even in very remote areas. Currently, three methods are commonly employed to derive glacier-surface velocities: interferometry of synthetic aperture radar (SAR) imagery, SAR tracking techniques, and cross correlation of optical satellite images. 3807 D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 Velocity measurements by interferometry of SAR imagery (InSAR) may achieve high accuracies, but require that coherence between the images is not lost due to modification of the glacier surface by, e.g., melting or snowfall (Strozzi et al., 2002; Trouvé et al., 2007). This requirement, together with limitations regarding the resolvable displacement gradients, result in InSAR-derived velocity measurements that are typically constrained to time spans of 1, 3 or 6 days (e.g., Massonet & Feigl, 1998; Joughin et al., 1996). Thus, the obtained velocity data may be representative only for the observation period and an extrapolation to annual velocities is difficult. Offset tracking in SAR imagery (Michel & Rignot, 1999; Joughin, 2002; Strozzi et al., 2002) is similar to cross correlation of optical satellite imagery (Lucchita & Ferguson, 1986; Bindschadler & Scambos, 1991). The basic approach is to track features from one scene to another and to calculate their velocity given the temporal separation and the measured displacement. In the case of SAR images, this can be done using either the intensity or coherence of the complex radar images (Strozzi et al., 2002). Compared to InSAR, tracking techniques using SAR images are more useful for measuring flow velocities over longer periods. However, a general drawback of SAR imagery in steep mountainous terrain is the high incidence angle of the sensor, which may inhibit visibility of the target glacier, and require very accurate DEMs to correctly orthorectify the measurements (Trouvé et al., 2007). Using optical satellite imagery, the detail and accuracy of the measurements is largely limited by the ground resolution of the sensor, and by the ability to precisely co-register images acquired at different dates. The latter task is usually the most difficult and has led to inaccuracies on the order of 1 pixel, i.e.,15 m if ASTER imagery were used (Kääb, 2005; Stearns & Hamilton, 2005). Further errors may arise from changes in the satellite attitude during scanning of the images (Van Puymbroeck et al., 2000), and from an inaccurate DEM during orthorectification using a rigorous model (e.g., Toutin, 2004). A principle drawback of optical imagery is the dependency on cloud-free conditions. In summary, velocity measurements by InSAR are most appropriate for analyzing very short time scales, i.e., days, or where extrapolation to longer time scales is justified, e.g., in ice sheet studies (Joughin et al., 2002). Feature tracking, using SAR or optical imagery is more appropriate for analyses over longer periods. Although limited by cloud cover during image acquisition, cross correlation of optical imagery provides a quick and efficient way of measuring glaciersurface velocities. Importantly, a huge and global archive of optical images from glaciers already exists and new images are continuously acquired. In order to achieve the measurement accuracy required to infer, e.g., annual velocity variations, the co-registration requires high accuracy and errors due to attitude effects or inaccurate DEMs need to be minimized. Here, we evaluate the potential and the limits of a new application for orthorectification, co-registration and correlation of optical imagery, COSI-Corr (Co-registration of Optically Sensed Images and Correlation; Leprince et al., 2007), to measure glacier-surface velocities in mountainous terrain. We provide guidelines to improve the accuracy of the measurements and to assess their quality without available ground-truth data. This includes correction of offsets in the displacement maps due to attitude effects and due to elevation errors in the DEM. The methodological principles are applicable to a wide variety of optical satellite imagery and are demonstrated here using ASTER images. We have studied the glaciers in two Himalayan regions: Khumbu in Nepal and Garhwal in India, where the glacier shrinking is observed. First, we demonstrate the methodological principles, including quality assessment, on the relatively slow Khumbu glacier at Mount Everest. Second, we investigate and model displacement errors induced by systematic elevation errors in the SRTM-based DEM, at the Gangotri glacier group in Garhwal. In a further step, the recent velocity history of Gangotri glacier, situated in the headwaters of the Ganges, is analyzed to demonstrate the capabilities and the limits of the method to monitor glacier dynamics 2. Methods and data Table 1 presents the imagery analyzed in this study, along with details on the acquisition parameters. Although we generally avoided Table 1 List of the ASTER scenes used in this study Region Granule ID Date [yyyy-mm-dd] Sun azimuth [degree] Sun angle [degree] Incidence angle [degree] Orientation [degree] Cloud covera [%] Khumbu (case study 1) ASTL1A 0009280513510312080 ASTL1A 0010140513270106251 ASTL1A 0112200502290201111 ASTL1A 0210040500380210261 ASTL1A 0211210500340212070 ASTL1A 0301080500160303170 ASTL1A 0310230459290311050 ASTL1A 0410090458390410220 ASTL1A 0410250458240411040 ASTL1A 0411100458190411210 ASTL1A 0511130458410511190 ASTL1A 0511290458400512020⁎ ASTL1A 0512060504390512090 ASTL1A 0512150458320512180 ASTL1A 0602010458090602040 ASTL1A 0701190459340701220 ASTL1A 0109090542130109210 ASTL1A 0310100529250310220 ASTL1A 0310100529340310220 ASTL1A 0407240529140408100 ASTL1A 0508190534580508220 ASTL1A 0510150528360510180 ASTL1A 0609230535100609260 ASTL1A 0610090534580610120⁎ ASTL1A 0611100535050611130 2000-09-28 2000-10-14 2001-12-20 2002-10-04 2002-11-21 2003-01-08 2003-10-23 2004-10-09 2004-10-25 2004-11-10 2005-11-13 2005-11-29 2005-12-06 2005-12-15 2006-02-01 2007-01-19 2001-09-09 2003-10-10 2003-10-10 2004-07-24 2005-08-19 2005-10-15 2006-09-23 2006-10-09 2006-11-10 155.78 161.76 160.96 152.76 162.48 157.48 158.65 154.41 158.11 160.38 161.12 161.18 162.41 160.27 151.87 154.56 149.10 156.13 155.70 116.65 133.17 157.07 151.63 158.14 163.20 57.51 52.29 36.18 54.87 40.26 36.44 48.60 52.87 47.51 42.70 41.93 38.58 37.35 36.29 39.99 37.74 60.91 49.64 50.21 68.37 65.31 47.74 55.82 50.61 40.39 −2.870 0.022 0.025 −2.829 −0.041 −0.030 0.019 0.022 −2.873 −1.480 0.022 −0.019 8.588 0.016 −2.876 −2.867 5.699 −5.727 −5.727 −8.586 5.729 −8.583 2.878 5.729 2.873 9.26 9.26 9.26 9.26 9.26 9.26 9.26 9.26 9.26 9.26 9.26 9.26 9.31 9.26 9.26 9.26 9.56 9.56 9.51 9.56 9.56 9.56 9.56 9.56 9.56 63 70 43 49 36 48 25 72 77 55 47 45 76 43 40 67 52 44 13 40 87 69 52 62 57 Garhwal (case study 2) All given data were extracted from the metadata of the images. The orientation measures the angle between the along-track direction and North in a clockwise direction. The images that were used as the master images in the co-registration procedure are marked with an asterisk (⁎). a The listed cloud cover is taken from the images metadata and usually overestimates the true cloud cover. 3808 D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 scenes with heavy cloud and snow cover, we included a number of less than optimal scenes to test their suitability for velocity measurements. The different steps of our approach are organized in two work flows and presented in Fig. 1. The first group of tasks comprises orthorectification, co-registration, and correlation of the satellite imagery, followed by post-processing of the correlation results using COSI-Corr. COSI-Corr is a new software package that has originally been developed for the detection of coseismic displacement (Leprince et al., 2007, 2008; available for download from the Caltech Tectonics Observatory website, http://www.tectonics.caltech.edu). The software package is an IDL-based module for the remote-sensing platform ENVI© by RSI. The application allows processing of aerial as well as satellite imagery from the SPOT, ASTER, and Quickbird sensors. A detailed description of the methodological background and COSI-Corr can be found in Leprince et al. (2007), and applications in Leprince et al. (2008) and Avouac et al. (2006). The second group of tasks is related to data filtering and assessing the quality of the results. In case of more than one correlation, i.e., more than two ortho-images, further steps may involve the comparison and the combination of the acquired data. 2.1. Orthorectification, co-registration and sub-pixel correlation of satellite images using COSI-Corr The orthorectification procedure relies on the automatic generation of ground control points (GCPs). A precise set of GCPs is generated from a raw image (slave), with respect to an already orthorectified image (master), by iteratively refining an initial rough selection of manually defined tiepoints. Image patches from the raw slave image are orthorectified and their misregistrations with the master image are estimated from correlation. A precise set of GCPs is produced when the misregistration measured at each patch converges to a minimum. Importantly, generating GCPs is independent of any ground data by using a shaded image of the DEM as the first orthorectified master. The first orthorectified image produced will then become the new master for subsequent slave images. This approach is globally applicable, wherever DEMs are available. However, the DEM needs to be free of voids, which is a common problem in mountainous terrain. Smaller gaps can be safely interpolated using standard methods while larger patches should be replaced with other data sources, as described in numerous studies (e.g., Luedeling et al., 2007; Crippen et al., 2007). Alternatively, SRTM tiles from many mountainous regions in the world, where most of the largest voids have been patched with data from topographic maps, are publicly available from Jonathan de Ferranti (http://www.viewfinderpanoramas.org). Such DEMs have been used in this study. Once a set of precise GCPs has been produced, the mapping matrices that associate ground coordinates with raw pixel coordinates are computed. They define the resampling grid from the raw image to the orthorectified image (Fig. 1A). Special care is brought to the resampling operation in order to avoid the introduction of aliasing in the orthorectified image. Horizontal ground displacements are retrieved from the sub-pixel correlation of multitemporal orthorectified images (Fig. 1B). Image correlation is achieved with an iterative, unbiased processor that estimates the phase plane in the Fourier domain. This process leads to two correlation images, each representing one of the horizontal ground displacement component (East–West and North–South), and to a Signalto-Noise Ratio (SNR) for each measurement, assessing the confidence of the results. In a typical process, images are wrapped onto the topography within the DEM resolution, and co-registered in pairs with 1/50–1/20 pixel accuracy, allowing for the measurement of horizontal offsets with an accuracy on the order of 1/20–1/10 of the pixel size. All data produced for this study have been obtained using ASTER band 3N 15 m resolution images. To allow the measurement of large displacements without losing resolution on the displacement fields, the COSI-Corr multiscale correlation analysis was performed using a window size of 128 down to 32 pixels. Steps of 4 pixels between adjacent correlations yielded ice-flow velocity maps sampled at every 60 m. Fig. 1. Processing chain of the applied method to derive accurate glacier-surface velocities. The first work flow comprises the orthorectification and co-registration of multitemporal satellite images (A), their correlation (B), and post-processing (C) to improve the accuracy of the displacement measurements. These steps were done using ENVI© with COSI-Corr. The correlation results are filtered (D) and then checked for their consistency using streamlines (E), stacked profiles (F), and strain maps (G) in the second work flow. D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 2.2. Post-processing procedures 2.2.1. Removal of residual attitude effects Data on the roll, pitch, and yaw of the satellite during image acquisition come with the imagery's metadata, and are accounted for during orthorectification. However, the ASTER sensor samples the attitude information not frequently enough to allow for full compensation of the resulting image distortions (Teshima & Iwasaki, 2008). As a result, the correlation maps of two ortho-images will show wave artifacts in the across-track direction of the image (cf. Fig. 4). A gentle long wavelength distortion in the along-track direction is attributed to focal plane distortions, e.g., spherical aberration from the optical system or distortion of the CCD sensor (Leprince et al., in press). Such systematic distortions can be removed using postprocessing tools within COSI-Corr. The possibility to remove these artifacts depends on the fraction of visible, stable ground, i.e., ground that does not involve any glacier flow, in the two ortho-images. Generally, the higher the amount of stable and visible ground, the better the possibilities of removing attitude effects. However, distortions resulting from attitude effects may be obscured when other distortions are present, e.g., due to inaccurate DEMs. 2.2.2. Removal of DEM-related errors Although COSI-Corr was explicitly designed for correlating satellite images irrespective of their incidence angles, different incidence angles may lead to distortions in the orthorectification in case of an inaccurate DEM. As these distortions cannot yet be corrected a priori, i.e., during orthorectification, they will be transferred to the displacement maps. In our case studies, errors were most prominent in the E–W displacement maps, as the ASTER sensor can only be inclined in the across-track direction and the orbital path of the carrying satellite TERRA is only a few degrees off north over sub-tropical latitudes. If we assume that all pixels in an image have a comparable sight angle that is well approximated by the instrument incidence angle, the measured ground disparity D relates to the incidence angles of the correlated scenes, θ1 and θ2, and to the elevation error of the DEM, h, by D ¼ h⁎ðtanðθ1 Þ − tanðθ2 ÞÞ: ð1Þ The disparity increases with the difference in incidence angles and the elevation error of the DEM (see Fig. 2). As the SRTM data is the principal source of DEMs in many studies, it is useful to assess any 3809 systematic errors that can be modeled to improve the accuracy of the displacement measurements. It has been shown in earlier studies that elevation errors of SRTM-based DEMs contain a component which linearly increases with terrain slope, and another which depends on terrain aspect (Bourgine & Baghdadi, 2005; Gorokhovich & Voustianiouk, 2006; cf. Toutin, 2002a). The dependency on terrain aspect is presumably related to the orbital path of the Space Shuttle and to the look direction of the antenna (Bourgine & Baghdadi, 2005). Accordingly, elevations of foreslopes, i.e., with a northwesterly aspect, are generally underestimated and elevations of backslopes, i.e., with a southeasterly aspect, are generally overestimated. Because the orthoimages and correlation maps are well co-registered with the DEM used for orthorectification, the ground disparities can be compared with the topography to produce a model for correcting the displacement errors. We found that the residual displacement error, ε, can be estimated with the model e ¼ K ⁎ s ⁎ cosða þ uÞ þ z ð2Þ where s is the slope of the topography surface, a is the topography aspect, and K, φ, and z are constants to be determined from, e.g., a least squares procedure. In all cases we investigated, φ was around 1.3 rad, i.e., 75°, which implies that the largest offsets occur at aspects of approximately 105 and 285° (see Table 2). K can be interpreted as the maximum offset among all aspects, per slope radian. In this study, the absolute value of K, for the E–W displacement, was always around 13 m/rad, i.e., about 23 cm per degree slope angle. The last term, z, is not related to the DEM-error but may be regarded as the mean error due to attitude effects. This term could be set to zero if the correlation results, after correcting for DEM-error effects, allow removal of the attitude effects with the destriping tool (cf. Table 2). In some cases (see below), this is not possible due to residual noise in the correlation map, which stems from (1) inaccurate slope and aspect values and (2) erroneous sampling of miscorrelations or moving ground for estimating the parameters K and φ. Before fitting Eq. (2) to the displacement, aspect and slope data, we used a signal-to-noise ratio threshold of 0.99 and a data range between −20 m and +20 m for the E–W and −10 m and +10 m for the N–S displacement to minimize noise and erroneous sampling. As an alternative to using SRTM-based DEMs, one could use the ASTER images from the 3N and 3B bands to construct the DEM used for orthorectification. However, it should be noted that the above mentioned attitude effects, as well as steep slopes, shadows, clouds, and snow fields will cause problems in the DEM generation. Consequently, ASTER-derived DEMs from steep, mountainous terrain are usually associated with many gaps and large errors (e.g., Toutin, 2002b; Kääb, 2002; Eckert et al., 2005; Fujisada et al., 2005), and thus we preferred to use the SRTM-based DEM. 2.3. Data filtering Fig. 2. Effect of DEM-error on displacement measurements. Assume a pixel p1 from an image I1 acquired at a date t1 sees the ground point M, and a pixel p2 from an image I2 acquired at a date t2 sees the same point M on the ground, and that both images are orthorectified and co-registered according to a DEM with an elevation error h. For simplicity, it is assumed that locally, around the ground point M, the topography and the elevation error are well approximated by constants. θ1 and θ2 are the angles between the line of sight of pixels p1 and p2, and the vertical. When the orthorectified images I1 and I2 are correlated, a disparity D = δ1 − δ2, induced by the elevation error h, is measured. Once all systematic errors have been removed, the measurements were filtered to exclude miscorrelations and to identify reasonable correlations obscured between miscorrelated patches (Fig. 1D). Excluding measurements with a low signal-to-noise ratio is a starting point to quickly filter the displacement maps. However, this does not exclude all miscorrelated points, and we have found that in addition, a simple directional filter is very efficient in getting rid of most remaining miscorrelations (e.g., Kääb, 2005). This was done by defining the flow direction from flow features on the glacier surface in the ortho-images and allowing for some deviation, e.g., of up to 20°. An additional filter is applied to the magnitude of the displacement to acknowledge that velocities do not change abruptly, but rather gradually. Both filter procedures need to be applied with variable parameters (e.g., directions, sizes, and thresholds) on different patches of the glacier and thus require some manual tuning. Overlaying the displacement field in form of vector arrows on one of the ortho- 3810 D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 Table 2 Details on the error evolution during post-processing of the correlations used in the study of the recent velocity history of the Gangotri glacier Correlation details Residual offset [m] Ortho 1 Ortho 2 Time span [a] Inc. angle diff. [degree] Aug-05 Sep-06 1.08 2.85 Aug-05 Oct-06 1.17 0.00 Aug-05 Nov-06 1.25 2.86 Jul-04 Oct-05 1.25 0.00 Jul-04 Aug-05 1.08 14.32 Oct-03 Jul-04 0.75 2.86 Oct-03 Aug-05 1.84 11.46 Oct-05 Oct-06 1.00 14.31 Oct-05 Sep-06 0.92 11.46 Oct-05 Nov-06 1.08 11.46 Oct-03 Oct-06 3.00 11.46 Sep-01 Oct-03 2.08 11.43 Sep-01 Aug-05 3.95 0.03 Raw E–W N–S E–W N–S E–W N–S E–W N–S E–W N–S E–W N–S E–W N–S E–W N–S E–W N–S E–W N–S E–W N–S E–W N–S E–W N–S Parameters for the correction model DEM-error corr. Attitude-effect corr. Mean Std Mean Std Mean Std −1.140 1.723 −0.741 −0.197 −1.021 −1.065 0.200 −0.062 0.132 1.465 −1.406 1.337 1.593 0.175 0.086 0.123 −0.173 −0.742 0.090 0.208 1.410 0.070 −1.560 −0.030 1.052 0.111 3.756 4.239 3.921 4.336 4.110 4.428 3.705 3.978 5.251 4.157 4.189 3.558 5.182 4.573 5.252 3.576 5.315 3.796 5.367 3.733 5.092 3.668 5.066 4.115 4.023 4.232 \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ 4.446 3.487 \ \ 5.033 4.536 4.546 3.484 4.972 3.768 4.844 3.485 4.449 3.275 3.612 2.851 \ \ −0.255 −0.752 −0.070 −0.038 −0.104 −0.221 0.039 −0.016 −0.074 0.191 0.116 0.308 0.230 0.385 −0.036 \ 0.045 −0.021 −0.008 0.226 \ \ 0.003 −0.007 0.299 −0.218 3.394 3.720 3.499 3.732 3.682 3.989 3.489 3.600 4.070 3.179 3.366 3.276 4.842 4.251 4.371 \ 4.614 3.591 4.428 3.423 \ \ 3.487 2.811 3.268 3.269 0.114 0.109 \ \ 0.174 0.145 −0.007 0.049 0.066 −0.091 0.009 0.113 0.135 0.107 −0.060 −0.033 \ \ K [m/degree] Phi [degree] z [m] \ \ \ \ \ \ \ \ −0.239272 0.048 \ \ −0.214183 0.043385 −0.243768 0.046721 −0.238114 0.049008 −0.232179 0.055981 −0.223690 0.036839 0.253199 −0.044595 \ \ \ \ \ \ \ \ \ \ 63.58685 68.48564 \ \ 78.20873 86.78591 71.33324 86.48574 70.74831 78.63152 73.25666 59.63550 72.83800 76.66753 69.16282 76.61620 \ \ \ \ \ \ \ \ \ \ −1.068 0.079 \ \ −4.585 −0.040 −1.020 0.002 0.240 0.955 −0.321 −0.164 −3.540 −0.040 3.971 −0.100 \ \ When the differences in incidence angles were low, corrections of DEM-induced errors were not necessary. Residual offsets were determined from all displacement data in a range between −10 m and +10 m. Thus, slow moving glacier ice has also unwillingly been sampled, and the residual offset estimates should be regarded as upper bounds. images helped to identify whether the results were consistent with the flow features on the glacier surfaces. We designed an interface in MATLAB© that allows for a quick definition of thresholds and patch sizes to apply the filters. glacier surface. The reference frame is the local flow direction. With this suite of tests, we determined whether the correlation procedure was stable and we produced consistent results that are supported by flow features on the glacier surface. 2.4. Quality assessment and validation techniques 3. Study area The lack of ground-truth velocity measurements hampers simple evaluation of remotely sensed measurements in most cases. Yet, in order to assess the quality of the measurements, we designed a number of tests to check the consistency of the results with regard to the displacement direction, magnitude, and gradient. These include (1) a test of the displacement direction by using the displacement field to construct streamlines, i.e., displacement paths, which can be checked against flow features on the glacier surface in the orthoimages (Fig. 1E); (2) a test of the displacement magnitude by comparing the sum of incremental measurements (e.g., the sum of displacements measured from images between 2001–2002, 2002– 2003 and 2003–2004) with a displacement measurement over the complete observation period (i.e., 2001–2004) (Fig. 1F); (3) a check of the displacement gradients by overlying the ortho-images with strainrate maps calculated from the displacement data (Fig. 1G), using the method by Nye (1959) as shown in studies by Bindschadler et al. (1996). For the calculation of strain rates, only filtered displacement values have been used and small gaps in the displacement maps have been linearly interpolated. Furthermore, in order to suppress smallscale dynamics and noise in the strain rates, the displacement maps have been smoothed with a 5 × 5 pixel convolution filter (Bindschadler et al., 1996). An error estimation of the strain-rate calculations was performed by bootstrapping (n = 1000) the calculations using the E–W and N–S displacements with added uncertainties. The uncertainties have been randomly drawn from a normal distribution described by the residual error over stable ground. The resulting strain-rate maps describe the longitudinal, transverse, and shear-strain rates over the Currently, approximately 116,000 km2 of mountainous terrain are glacierized in South and Central Asia (Dyurgerov & Meier, 2005), making this region the largest glacierized continental area outside the polar regions. Despite the great number of glaciers in the Himalaya and Karakoram and their important role for water supply to the region, glaciological data are surprisingly limited (e.g., Wagnon et al., 2007). The available measurements of glacier areas and mass-balance calculations have shown that glaciers in the Asian highlands are generally retreating (Mayewski & Jeschke, 1979; Dyurgerov & Meier, 2005), in some cases at high rates, like the Parbati glacier in India, retreating at almost 52 m/year (Kulkarni et al., 2005). Conversely, some glacier advances have been observed in the eastern Himalaya and the Karakoram, which have been linked to increased precipitation (Liu et al., 2006) and/or decreased summer temperatures (Hewitt, 2005; Fowler & Archer, 2006). Because of the low latitudinal position between 27 and 37° N, Himalayan glaciers usually occur at elevations of more than 4 km, although some descend to elevations of less than 3 km. Such advances to relatively low altitudes are commonly thought to be driven by a high amount of supraglacial debris cover that shields the ice from ablation, lowering the accumulation-area ratios compared to that of debris-free glaciers (Benn et al., 2003). The debris cover is an important feature for deriving surface velocities from optical satellite imagery as it creates and preserves pronounced surface morphology over relatively long timescales (Luckman et al., 2007). However, the correlation procedure tends to fail when illumination conditions are grossly different between scenes. During summer, frequent cloud D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 cover due to the Indian Southwest Monsoon limits the choice of suitable satellite scenes. In our study we have chosen the Mount Everest region, Khumbu, in the Nepalese Himalaya, and the Gangotri glacier group, Garhwal, in the Indian Himalaya. We selected these sites because they hold abundant glaciers of different sizes that are important water resources (e.g., Singh et al., 2006) and some of them, due to recent downwasting, are prone to catastrophic outburst flooding (e.g., Cenderelli & Wohl, 2001; Kattelmann, 2003), making them prime targets for monitoring strategies. The high elevation sectors in both regions are characterized by a moderately wet monsoonal climate, with more influence of the Winter Westerlies in Garhwal. 4. Results 4.1. Case study 1: Khumbu Himal, Nepal Fig. 3 shows an ASTER ortho-image from the Mount Everest region, acquired in November 2004, and a displacement map produced by correlation with another ortho-image from November 2005. The acquisition setting of both ASTER scenes with identical near-vertical incidence angles, similar shading, absence of clouds, and only limited snow cover, provide ideal conditions for orthorectification and correlation (see Table 1). 4.1.1. Accuracy improvement Well-identifiable stripes in the E–W displacement map are due to attitude variations and are a first sign of low noise and successful orthorectification (Fig. 4). The stripes have been removed with the COSI-Corr destriping tool. This has improved the accuracy of the measurements as is shown in Fig. 4. Before the correction, the residual displacement in the E–W direction, as measured from all data points lying within −10 m to +10 m, had a mean value of −0.63 m and a standard deviation (1σ) of 3.16 m. After removing the distortions in the line direction of the image, the residual displacement decreased to a mean value of −0.11 +/− 2.52 m. Further removal of the more gentle distortion in the column direction improved the accuracy only somewhat to a mean of −0.05 +/− 2.41 m. Most likely, optimal results 3811 from the destriping procedure would be achieved if the destriping model was defined using stable ground only. However, this would have been a laborious task, and we found that simply thresholding the displacement map to a range that encompasses the undulations due to attitude effects, e.g., −10 m to +10 m, works well enough to remove any systematic undulations. Indeed, most of the glaciers have moved by more than 10 m during our study period, and thus most of the iceflow related measurements are discarded from simple thresholding. 4.1.2. Filtering After removal of obvious systematic distortions in the displacement images, the displacement measurements over the glacier area have been filtered to eliminate miscorrelations. This approach is used on Khumbu glacier (N27.9806, E86.8766), which is an intermediatesize glacier (16.5 km length) located southwest of Mount Everest. Based on an analysis by Luckman et al. (2007), the glacier appears to be stagnant over its lowermost 2 km. As was already apparent in the displacement map in Fig. 3B, the correlation procedure failed in certain parts and returned erroneous displacements (see Fig. 5). This is particularly the case in the steep portions of the glacier where the velocity gradient, and thus deformation of the glacier surface, is large. Another problem that may not be apparent at first sight is artificial displacement due to moving shadows (e.g., Berthier et al., 2005). If the sun angles are different in the scenes to be correlated, the correlation procedure will possibly detect the shifting shadows and record an artificial displacement. In order to exclude such miscorrelations, we have filtered the data over the area of the glacier as described in Section 2.3. The result of the filter procedure in the central part of Khumbu glacier is shown in Fig. 5. Most of the obvious erroneous vectors have been discarded using the directional filter (black arrows). The magnitude filter discarded another group of displacement vectors that were pointing in the correct direction, but showed anomalously high or low displacement values (red arrows). In this case, we applied the filters on patches of up to 1 km2, depending on changes in flow direction and magnitude, and allowed for +/−20° deviation from the defined flow direction. The magnitude filter was applied more variably according to nearby, well identifiable velocities, usually within a range of +/−30 m/a. Clearly, Fig. 3. Ortho-image (A) acquired on Nov. 10th 2004 and displacement map (B) from the Mount Everest region, Nepal. The displacement map was produced by correlating the orthoimage in A with another ortho-image acquired on Nov. 29th 2005. The displacement values were normalized to annual velocities. The subsets in A are shown in Figs. 5 and 7 and the velocity along the profile in B is displayed in Fig. 6 (short profile) and Fig. 8 (long profile). 3812 D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 Fig. 4. Correction of attitude effects (A) and sensor distortions (B) in the E–W displacement component of the correlation shown in Fig. 3B. The corrections are calculated from the mean residual offset in the column- (A) and in the line direction (B) of the image. The scatter plots depict the individual offsets of a randomly sampled set of 10,000 pixels, and the histograms show the entire population arranged in 0.02 m offset bins. Fig. 5. Velocity field of the central part of Khumbu glacier derived from the correlation shown in Fig. 3B. The ortho-image in the background is from Nov. 10th 2005. The arrows depict displacements of more than 1.5 m/a. Through filtering the data by direction, most miscorrelations are discarded. Further filtering by magnitude removed measurements pointing in the right direction, but that showed anomalously high or low velocities (red arrows). Green arrows mark the filtered velocity vectors that are consistent with flow features on the glacier surface. Streamlines are shown in white and were constructed using the retrieved velocity vectors. D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 unless filtering is performed very carefully with tight thresholds and on small patches, erroneous results may survive and correct results may be discarded, e.g., some displacement vectors at the edges of Khumbu glacier in Fig. 5. However, the proportion of these cases among the entire population of retrieved data is usually very small, even if the filtered patches are relatively large. 4.1.3. Quality assessment The first consistency test using streamlines is applied on the lower part of Khumbu glacier, which has good and continuous data coverage (Fig. 5). The streamlines agree quite well with the flow features on the glacier surface. In the lower part of Khumbu glacier the streamlines are narrow, due to the confluence with a tributary glacier. A minor mismatch of the streamlines coming from the tributary glacier and the banding on the glacier surface (hardly visible in the figure) does not appear to be an artifact, as it is observable in all other correlations we performed. Instead, the mismatch appears to reflect a relative increase in ice discharge from this tributary compared to the upper Khumbu glacier. For the second consistency test we used four ortho-images from the years 2001, 2002, 2004, and 2005 (see Table 1). An example from the lower part of Khumbu glacier is given in Fig. 6. For this profile the raw, unfiltered data have been used to show the good agreement over most parts of the profile. The displacements derived between 2002 and 2004, and all time spans encompassing this period, show some suspicious velocity variations in the center of the profile. From visual inspection it was found that these variations are due to the enlargement of a supraglacial pond, where the retreat of the bounding margins caused an apparent reduction in velocity at the up-glacier side and an apparent increase in velocity at the down-glacier side of the pond. The stacked profiles show that the magnitudes of the displacement measurements agree very well in the upper part of the profile but contain larger scatter in the lower part, where surface degradation due to melting is higher. Furthermore, at lower displacement values, the distorting influence of the noise increases, especially regarding the displacement direction. These poor quality data points have been excluded using the filter procedure. It should be noted that the measured displacements are always straight and may thus lead to 3813 an underestimation of the true displacements if their paths were curved. However, as the magnitudes of the displacement vectors are generally small compared to the local curvatures of the flow, the displacement paths are well approximated linearly. Problems may occur when measuring flow in strong bends over longer time spans. For the third consistency test, we calculated strain rates from the displacement field and examined its compatibility with regard to the glacier surface, e.g., the occurrence of crevasses. However, the correlation failed in the central part of Khumbu glacier where crevasses are formed. Nevertheless, an examination of the pattern of strain rates still allows identification of unexpected displacement gradients. Fig. 7 shows the components of the calculated surface strain rates over Khumbu glacier and the error in longitudinal-strain rates. While most of the glacier is characterized by moderately low strain rates, some areas stand out with much higher strain rates. First, in the highest parts on Khumbu glacier where the velocity data were retrieved, the glacier slows down considerably, which causes high values of negative strain rates, i.e., shortening rates. This happens just below a steep part along the glacier profile, where numerous crevasses have formed, and presumably closed again. Second, approximately 500 m west, along-flow shortening reaches a peak at the confluence with a tributary glacier coming from the north. When looking at the velocity vectors and streamlines in Fig. 5 and at the transverse strain rates in Fig. 7, it is apparent that ice near the edge of Khumbu glacier has divergent flow towards the tributary glacier. Consequently, the tributary ice, which flows with velocities of less than 3 m/a near the confluence, is being pushed aside and not incorporated into the main ice stream of Khumbu glacier. Therefore, the contribution of ice from the tributary glacier appears to be reduced, which causes the Khumbu glacier to expand laterally. Newly formed crevasses with a NW–SE orientation can be seen in high-resolution satellite images (e.g., in Google Earth©), supporting this conclusion. In the upper part of the covered area of Khumbu glacier, shearstrain rates at the glacier margins are high and of opposite sign, as would be expected. In the lower part, where surface velocities as well as velocity gradients across the width of the glacier are low, shearstrain rates are lower too. The error on the longitudinal strain rates Fig. 6. Stacked displacement profiles from the lower part of Khumbu glacier (see Fig. 3B for footprint). The data points depict the displacement over the time span given in the figure's legend. The red line (blue line) is the sum of the 2 (3) displacements denoted by ‘+’ and ‘x’ (and ‘Δ’). The red (blue) circles depict the displacement measured over the same time periods as covered by the red (blue) lines. The difference in displacement and direction of displacement between the open circles and the lines is shown in B and C, respectively. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.) 3814 D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 Fig. 7. Longitudinal (A), transverse (B), and shear (C) strain-rate maps and the error on the longitudinal strain rates (D) over Khumbu glacier derived from the filtered velocity vectors shown in Fig. 5. The distribution and magnitude of transverse and shear-strain rate errors look similar to the longitudinal strain-rate errors and are not shown. See text for details on the calculation of the strain rates and strain-rate errors. (Fig. 7D) is highest in the regions of low velocities, as the flow direction is strongly affected by the noise, which consequently results in considerable scatter of the velocity gradients. The errors on the transverse and shear-strain rates are similar to that of the longitudinal-strain rates, and are therefore not shown. 4.1.4. Data combination: continuous velocity profile from the Khumbu glacier Because most correlation maps contain some areas where the procedure failed or returned erroneous data, it may be useful to combine the results from several correlations to enhance the spatial coverage across a glacier. A comparison of the filtered velocity measurements (not shown) has yielded very similar results throughout the observation period, from 2001 to 2007. In order to arrive at a continuous velocity profile of Khumbu glacier, we extracted the displacement data of all correlations we performed along a profile that extends from the highest point in the accumulation zone, down to the toe of the glacier (Fig. 3B). We applied our filtering procedure, with the same parameters, simultaneously to all displacement maps to only extract the meaningful data from our profile (Fig. 8). A large data set of 22 displacement maps enabled us in this case to produce a relatively well-constrained velocity profile, even though the results of the correlations are not equally good over the whole glacier. The lower part is especially consistent and the standard deviation among data points from different displacement maps is well below 5 m/a. In the central part, where the glacier flows over steep topography and attains high velocities and strain rates, strong surface modifications between the images complicated the correlation procedure and resulted in miscorrelations. It should be emphasized that the combination of velocity measurements from different time periods is only possible when the glacier shows no signs of velocity change over the period of observation. This condition has to be examined, e.g., using velocity profiles, before compiling the data. 4.2. Case study 2: Garhwal Himalaya, India The Gangotri glacier group is situated in western Garhwal, India, and forms part of the headwaters of the Ganges. The Gangotri glacier is, with more than 30 km length, one of the largest glaciers in the Indian Himalaya. We obtained 9 ASTER scenes, covering a period from September 2001 to November 2006, that are characterized by differences in incidence angle of up to 14° (see Table 1), which caused additional errors (see Section 2.2.2). Fig. 8. Continuous velocity profile of the Khumbu glacier derived from 22 correlations. While all data points (A) show considerable scatter in the central and upper part, the filtered data points (B) display less scatter as is shown in the plot of the standard deviation (C) of the data points from all the correlations. Note that not all correlations provide measurements over the entire glacier. 4.2.1. Correcting for DEM-related distortions Fig. 9 depicts the E–W displacement over Gangotri glacier and adjacent areas, derived from a correlation of ortho-images from October 2003 and October 2006 (see Table 1). The difference in incidence angles between the ortho-images is 11.5°. In Fig. 9A it is seen how displacement errors over stable ground produce an artificial shading effect, which highlights the dependence of the elevation error on terrain aspect. The variation of the mean E–W and N–S offsets with terrain aspect and slope angle is given in Fig. 10. We modeled the offsets using Eq. (2), with K = −12.817 m/rad, φ = 1.271 rad, and z = −3.54 m for the E–W component, and K = −2.111 m/rad, φ = 1.338 rad, and z = −0.04 m for the N–S component, determined from least squares adjustment. The applied correction improved the measurement accuracy to the degree that distortions due to the attitude effect became visible (Fig. 9). However, in this case we were not able to further correct the attitude effects more precisely as described in Section 2.2.2, due to a high level of noise and a small fraction of stable, correlated ground that could be used for defining the destriping D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 model. Thus, the negative z value in the E–W component represents the mean attitude effect which is biased towards higher values in the upper part of the image, where more stable, correlated ground exists. The correction improved the mean residual errors determined from all displacement values between −10 m and +10 m, from 1.41 +/− 5.1 (errors are 1σ) m to 0.13 +/− 4.4 m for the E–W component, and from 0.07 +/− 3.6 m to 0.11 +/− 3.2 m for the N–S component. Better results were obtained in cases where additional destriping was possible (Table 2). Nevertheless, given that distortions from DEM-errors increase linearly with slope angle, the impact on the derived glacier velocities is only small as glaciers mostly occur on low-gradient terrain. This is shown in Fig. 9C, which depicts the surface velocities of the Gangotri and the adjacent Chaturangi glaciers along profiles from different correlations. A measurement from September 2001 to August 2005, with no difference in incidence angle, is used as a reference. Although natural velocity variations may occur, these should be rather small due to the length of the observation period. The velocity difference between the uncorrected and corrected correlation (October 2003–October 2006) is small and almost not visible. Furthermore, the velocity measurements from the corrected correlation and the correlation from September 2001 to August 2005, yield very similar values. The results from the error modeling and removal of other correlations used in this study are given in Table 2. At incidence 3815 angle differences of more than 10°, DEM-induced errors were visible, modeled and removed. In most cases it was possible to correct the displacement maps for attitude effects after removal of the DEMinduced errors. The residual errors of the corrected displacement are similar to the residual errors of correlations with low incidence angle differences. Thus, the error removal was successful. 4.2.2. Data comparison: recent velocity history of the lower part of Gangotri glacier Velocity measurements from the correlations presented in Table 2 were used to investigate the recent velocity history of Gangotri glacier. For this purpose, we picked a profile along the central flow line of the glacier and plotted the annual velocity, with the associated errors given as shaded areas around some of the measurements, in Fig. 11. Over most of the profile, the annual velocity from October 2003 to July 2004 was faster than during the period from July 2004 to October 2005 (Fig. 11A). The difference is greater than the combined error on the measurements, and is therefore significant. The annual velocity from October 2003 to August 2005 rests in between the analyzed periods as would be expected. Whether or not this velocity difference is a true decrease in ice discharge over time or an effect of the sampled period, e.g., a seasonal effect, is not clear from this analysis. In order to elucidate the role of the seasonal coverage of the observation periods, we investigated annual velocities over different Fig. 9. Uncorrected (A) and corrected E–W displacement map (B) over the Gangotri glacier group derived from the correlation of ortho-images acquired on Oct 10th 2003 and Oct 9th 2006. The color-coding applies to pixels with displacement values in between +20 m and −20 m. Pixels with values outside this range are assigned the last color of the respective side of the color bar. The aspect dependence of the SRTM-based DEM-error produces an artificial hill shade effect in (A). After correction of DEM-induced offsets, distortions due to attitude artifacts become visible (B). Velocity profiles of the Gangotri and Chaturangi glaciers are shown in (C). For comparison, the uncorrected and corrected velocity measurements from the period 2003–2006 are plotted together with velocity measurements from the period 2001–2005, where topography-related artifacts are absent as the incidence angles during acquisition of the images were similar (see Table 2). 3816 D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 Fig. 10. Mean E–W (A) and N–S offset (B) from the correlation of the Oct 2003 and Oct 2006 ortho-images from the Gangotri glacier group (see Table 1) as a function of slope angle and aspect. The mean offsets were determined from the E–W and N–S displacement values ranging between −20 m and +20 m, and −10 m and +10 m, respectively. Slope angles were sampled in 1 degree intervals and aspect in 5 degree intervals. The resulting distribution was smoothed with a 5 × 5 convolution filter. At slope angles of more than 45° (not shown), fewer data points lead to more scatter. periods within one hydrological year, although only available from 2005 to 2006. Fig. 11B depicts a less obvious, but still visible difference in annual velocity when comparing a time period starting in August 2005 with one starting in October 2005. Velocities from October 2005 to October 2006 appear faster than velocities from August 2005 to October 2006. Very similar velocity profiles from two more correlations ending in late September 2006 lend additional credibility to the results. Importantly, the occurrence of the velocity difference spatially coincides with the larger velocity difference observed in the earlier time periods in Fig. 11A. The main difference of the seasonal coverage between the observation periods with slower and faster surface velocities is the extension of the slow velocity observations into the third quarter of the year, i.e., over late July to early October in the first case (Fig. 11A) and over late August to early October in the second case (Fig. 11B). Hence, the flow velocity during this time appears to have been relatively slower compared to the average velocity during the rest of the year. The larger velocity difference in the first case (Fig. 11A) can be explained if periods of slower velocities extended from July to October in both years, 2004, and 2005. Therefore, we conclude that the measured difference in annual velocities from 2003 to 2005 may be due to the same reason as for the observed difference in velocity during the period from 2005 to 2006. Several studies on alpine glaciers as well as outlet glaciers of ice sheets have shown that glacier flow velocities can vary significantly over daily to annual time scales (Bindschadler et al., 1977; Gudmundsson et al., 2000; Anderson et al., 2004; Bartholomaus et al., 2008). Such variations have commonly been attributed to melt water induced changes in subglacial hydrology that lead to variations in the speed of basal sliding. Importantly, many mountain glaciers show the highest flow velocities during spring to early summer, and before maximum ablation and proglacial stream discharge occurs (Willis, 1995; Fountain & Walder, 1998; Harper et al.; 2007). Such phenomena may explain the observed variations in flow velocity of the Gangotri glacier as well. During early summer, velocities may be higher as temperatures are high and melting occurs. However, as melting decreases from August to October, flow velocities may be slower than the annual average. Such inferences are supported by meteorological and hydrological studies near the terminus of Gangotri glacier (Singh et al., 2006, 2007). Therefore, we speculate that subsequent to peak melting and discharge in July/August, flow velocities decrease to slower than average levels. Hence, the observed decrease in average annual velocity from 2003 to 2005 may be the result of the observation period and may not reflect an overall decrease in flow velocity and ice discharge. 5. Discussion 5.1. Measurement errors All data in this study are reported as horizontal glacier-surface velocities or displacements. The data have not been converted to surface-parallel velocities. This can be easily achieved with the DEM Fig. 11. Recent velocity history of the lower part of the Gangotri glacier derived from the correlation of ortho-images from the years 2003 to 2006. The shaded rims around selected profiles give the one-sigma error range, calculated from the residual offsets (see Table 2). (A) shows the annual velocity during the period 2003 to 2006. Significant differences exist between the period 2003–2004 and 2004–2005 over the reach marked by a gray background. (B) shows the annual velocity from different observation periods from 2005 to 2006. Although of lower magnitude, a similar velocity difference is visible over the same reach as in A. The footprint of the profile is shown as Profile 1 in Fig. 9A. D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 used for orthorectification, as the topographic and kinematic data are well co-registered. However, such conversion does not account for the emergence velocity, which is the vertical velocity due to accumulation and ablation (e.g., Reeh et al., 2003). We estimated the residual errors on the measurements by analyzing the distribution of offsets with absolute values of less than 10 m. This means that “slow” moving ice is erroneously sampled, hence skewing the distribution to higher offsets. Applying an additional threshold of 0.99 to the signal-to-noise ratio map usually limits the data to low-relief areas. This results in much lower residual offset values on the order of zero mean and a standard deviation of 1 m. However, as we cannot assume that the residual errors on measurements of moving glaciers can be characterized by stable ground with high SNR values, we applied the rather conservative error estimation without using a high SNR threshold. Hence, the errors presented in this study should be regarded as upper bounds. As in many regions worldwide, the most accurate DEM available to date for high mountainous terrain is based on SRTM data, and thus problems associated with systematic topographic errors may occur frequently. Our method presents a means to model and correct for the resulting displacement errors. Yet, at slope angles greater than ∼ 45°, the model does not fit the offset data as good as at lower slope angles, due to large scatter and insufficient data points. Fortunately, glaciers occupy mostly low-gradient terrain, where such topography-related errors are small, providing good possibilities for correction. Furthermore, an advantage of ASTER imagery over most other sensors is that the incidence angle of the 3N band (VNIR), which should be used for velocity measurements, is always close to nadir, hence assuring small distortions due to topographic errors. Note that the presented error description and modeling is related to our use of an SRTM-based DEM specifically, as we observe the bias of the DEM, scaled by a function of the incidence-angle difference (Eq. (1)). The error modeling may be applied to all correlation results of ortho-images produced with SRTM-based DEMs, but the fitting parameters are specific to each correlation. 5.2. Comparison with SAR-derived velocity measurements When comparing our velocity measurements of Khumbu glacier with those obtained by Luckman et al. (2007) using InSAR and feature tracking, important differences emerge. First, the SAR-derived data show significant scatter over most of the profile, which is not seen in our data. Second, the measurements obtained from feature tracking and interferometry over two time periods each, differ considerably between the techniques and also between two periods using one technique. Problems with InSAR in the lower parts of the analyzed glaciers are acknowledged by these authors and more confidence is put on the data obtained from feature tracking. Luckman et al. (2007) also analyzed Kangshung glacier on the eastern slopes of Mt. Everest. In this area, they obtained results that are much more consistent with our data. However, the errors given by these authors are quite high, sometimes exceeding the velocity and on average reach 50% of it. These errors are thus too high to reliably detect velocity changes and assign them to natural causes and not to measurement problems. Despite these difficulties it should be noted that tracking techniques using SAR imagery may complement the use of optical imagery in areas where the glacier surface lacks visible features, e.g. in extensive accumulation troughs. 5.3. Other optical sensors Apart from ASTER, imagery from other satellites, as well as aerial photos, can also be used with COSI-Corr to measure ground displacement (e.g., Leprince et al., 2007). Satellite pour l'Observation de la Terre (SPOT) imagery in particular has proven useful in deriving glacier-surface velocities (Berthier et al., 2005; Leprince et al., 2008). Compared to ASTER, SPOT images have a more accurate attitude 3817 description (attitude variations are sampled at a higher rate), and do not usually require the correction of attitude effects in the displacement maps. However, as the incidence angles in SPOT images can be high, DEM-errors in steep terrain may cause larger distortions. It is also useful to know that the broader spectral bandwidth of the panchromatic sensor of SPOT (500–730 nm) has often led to high gains of many earlier SPOT images, when high mountainous terrain was not among the main target areas of satellite-data acquisitions. This resulted in saturated images that are useless for velocity measurements over snow-covered areas. More recent SPOT5 imagery is now adapted and provides images with lower gains over snowcovered mountains. It is not possible to process satellite images from the Landsat Thematic Mapper (TM) or Enhanced Thematic Mapper+ (ETM+) with sub-pixel accuracy. This is due to the unknown attitude variations of the satellite and the imaging system. Whereas SPOT and ASTER are pushbroom sensors, i.e., they scan across-track lines of 60 km instantaneously, TM and ETM+ sensors sample the ground along 16 across-track lines of 185 km (this is a whiskbroom system). Hence, attitude variations do not only occur in the along-track direction, but also in the across-track direction, which makes their removal virtually impossible. 5.4. Implications for glacier monitoring The suitability, global coverage, and low cost of ASTER scenes make this imagery a viable option among other alternatives to conduct large-scale and long-term monitoring campaigns of remote glacial systems (e.g., Kargel et al., 2005). In comparison with other sensors, the use of ASTER imagery, along with COSI-Corr, provides reliable results, as inherent problems with attitude effects and inaccurate DEMs can be solved. To successfully derive accurate glacier-surface velocities in mountainous terrain, a number of important points should be kept in mind. First, cloud cover should be low. However, when the master image has been successfully orthorectified, all other slave images require only three tie points to be accurately co-registered. Thus, even cloudy images with 3 patches (approx. 3 km × 3 km in size) of stable ground can be accurately co-registered. Importantly, thin, partly transparent clouds do not pose a problem. Therefore, even though cloud cover restricts the use of optical imagery to derive glacier-surface velocities, in many cases, images with even 50% of cloud cover or more can be used as long as the glacier of interest is visible. Second, images with grossly different snow-cover characteristics, such as winter and summer scenes, are difficult to correlate. The problem is not the snow cover itself, but the difference. That is why images from the same season usually work well, whether with or without snow cover. As the degree of snow cover is usually not identical between two images, parts on the glacier where the correlation procedure obtained poor results or failed are commonly encountered. Such data gaps may be filled with another correlation if images are present and the velocity did not change, as in the case of Khumbu glacier in the first case study. Third, surface modifications, like snow-cover change, complicate the correlation procedure. This applies directly to the resolvable time span and measurable velocities. When velocities are high, shorter time spans between the ortho-images lead to better results. For example, a surging glacier, which may flow at rates of several hundred meters per year, requires a narrow temporal separation of the images (in such a case, the use of InSAR may be more appropriate). When velocities are low, a longer temporal separation of the images is preferred, if surface degradation by melting or down-wasting does not interfere. Time spans exceeding one year also reduce the residual error when normalizing the results to annual velocities. For instance, we succeeded in measuring annual velocities of approx. 10 m/a on glaciers with little surface degradation in Garhwal. 3818 D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 Finally, we emphasize the need to properly understand the flow characteristics of the investigated glaciers before inferring any long-term trends. Our study of the velocity history of Gangotri glacier has shown that slight differences in the seasonal coverage of the observation periods may result in considerable velocity differences. Therefore, velocity analyses over different seasons should always be undertaken, if possible. 6. Conclusions In this study we have used the new application COSI-Corr to orthorectify, co-register, and correlate ASTER satellite imagery and derive surface velocities from glaciers in the Himalayas. We have shown how to minimize residual offsets of the displacement measurements due to attitude effects, and we have presented a way to detect, model, and correct for additional offsets induced by elevation errors of the SRTM-based DEM. Additionally, we developed techniques to check the quality and consistency of the results despite lack of ground control. The achieved measurement accuracies allowed us to detect seasonal velocity variations of 10–20 m/a in the lower part of the Gangotri glacier, in India. The results of individual correlations may be combined to enhance the spatial coverage across a glacier, which is particularly useful for synoptic studies aiming at continuous velocity profiles or maps from glaciers over large areas. We find that COSI-Corr is a powerful remote-sensing tool to perform detailed, as well as synoptic comparisons of glacier dynamics in remote, high mountainous terrains. Furthermore, the accurate coregistration of the ortho-images, the displacement maps, and the DEM used for orthorectification, provide the possibility to investigate links between surface features on the glacier, glacier dynamics, and topography. This may prove useful for future modeling studies that require tuning to recent conditions. Finally, our analysis of glacier dynamics using combined ASTER imagery and COSI-Corr has shown that this methodology is well suited to derive accurate, low-cost glacier-surface velocity measurements from remote regions where ground instrumentation is costly and difficult to implement. This is important in light of global warming and the need for water-management plans to take account of shrinking glaciers and the associated hazards. Acknowledgments This research was funded by a scholarship to D.S. as part of the German Research Council (DFG, Deutsche Forschungsgemeinschaft) graduate school GRK1364. ASTER imagery were provided to M.R.S. by the NASA Land Processes Distributed Active Archive Center User Services, Sioux Falls, SD, U.S.A. We greatly appreciate the constructive comments by three anonymous reviewers that helped improve an earlier version of this manuscript. References Aizen, V. B., Kuzmichenok, V. A., Surazakov, A. B., & Aizen, E. M. (2007). Glacier changes in the Tien Shan as determined from topographic and remotely sensed data. Global and Planetary Change, 56(3–4), 328−340. Anderson, R. S., Anderson, S. P., MacGregor, K. R., Waddington, E. D., O'Neel, S., Riihimaki, C. A., et al. (2004). Strong feedbacks between hydrology and sliding of a small alpine glacier. Journal of Geophysical Research, 109, F03005. doi:10.1029/2004JF000120 Avouac, J. P., Ayoub, F., Leprince, S., Konca, O., & Helmberger, D. V. (2006). The 2005, Mw 7.6 Kashmir earthquake: sub-pixel correlation of ASTER images and seismic waveform analysis. Earth and Planetary Science Letters, 249(3–4), 514−528. Bartholomaus, T. C., Anderson, R. S., & Anderson, S. P. (2008). Response of glacier basal motion to transient water storage. Nature Geoscience, 1, 33−37. Benn, D. I., Kirkbride, M. P., Owen, L. A., & Brazier, V. (2003). Glaciated valley landsystems. In D. J. A. Evans (Ed.), Glacial Landsystems (pp. 372−404). London: Hodder Arnold. Berthier, E., Arnaud, Y., Kumar, R., Ahmad, S., Wagnon, P., & Chevallier, P. (2007). Remote sensing estimates of glacier mass balances in the Himachal Pradesh (Western Himalaya, India). Remote Sensing of Environment, 108, 327−338. Berthier, E., Vadon, H., Baratoux, D., Arnaud, Y., Vincent, C., Feigl, K. L., et al. (2005). Surface motion of mountain glaciers derived from satellite optical imagery. Remote Sensing of Environment, 95, 14−28. Bindschadler, R. A., & Scambos, T. A. (1991). Satellite-image-derived velocity field of an Antarctic ice stream. Science, 252(5003), 242−252. Bindschadler, R. A., Vornberger, P. L., Blankenship, D. D., Scambos, T. A., & Jacobel, R. (1996). Surface velocity and mass balance of Ice Streams D and W, West Antarctica. Journal of Glaciology, 42(142), 461−475. Bindschadler, R., Harrison, W. D., Raymond, C. F., & Crosson, R. (1977). Geometry and dynamics of a surge-type glacier. Journal of Glaciology, 18(79), 181−194. Bourgine, B., & Baghdadi, N. (2005). Assessment of C-band SRTM DEM in a dense equatorial forest zone. C. R. Geoscience, 337(14), 1225−1234. Cenderelli, D. A., & Wohl, E. E. (2001). Peak discharge estimates of glacial-lake outburst floods and “normal” climatic floods in the Mount Everest region, Nepal. Geomorphology, 40, 57−90. Crippen, R. E., Hook, S. J., & Fielding, E. J. (2007). Nighttime ASTER thermal imagery as an elevation surrogate for filling SRTM DEM voids. Geophysical Research Letters, 34, L01302. doi:10.1029/2006GL028496 Dyurgerov, M. B., & Meier, M. F. (2005). Glaciers and the Changing Earth System: A 2004 Snapshot.: Institute of Arctic and Alpine Research, University of Colorado Occasional Paper No. 58, 117 pp. Eckert, S., Kellenberger, T., & Itten, K. (2005). Accuracy assessment of automatically derived digital elevation models from aster data in mountainous terrain. International Journal of Remote Sensing, 26(9), 1943−1957. Fountain, A. C., & Walder, J. S. (1998). Water flow through temperate glaciers. Reviews of Geophysics, 36(3), 299−328. Fowler, H. J., & Archer, D. R. (2006). Conflicting signals of climatic change in the Upper Indus Basin. Journal of Climate, 19(17), 4276−4293. Fujisada, H., Bailey, G. B., Kelly, G. G., Hara, S., & Abrams, M. J. (2005). ASTER DEM Performance. IEEE Transactions on Geoscience and Remote Sensing, 43(12), 2707−2714. Goldstein, R. M., Engelhardt, H., Kamb, B., & Frolich, R. M. (1993). Satellite radar interferometry for monitoring ice sheet motion: application to an Antarctic ice stream. Science, 262(5139), 1525−1530. Gorokhovich, Y., & Voustianiouk, A. (2006). Accuracy assessment of the processed SRTM-based elevation data by CGIAR using field data from USA and Thailand and its relation to the terrain characteristics. Remote Sensing of Environment, 104, 409−415. Gudmundsson, G. H., Bassi, A., Vonmoos, M., Bauder, A., Fischer, U. H., & Funk, M. (2000). High-resolution measurements of spatial and temporal variations in surface velocities of Unteraargletscher, Bernese Alps, Switzerland. Annals of Glaciology, 31, 63−68. Haeberli, W., Cihlar, J., & Barry, R. G. (2000). Glacier monitoring within the global climate observing system. Annals of Glaciology, 31, 241−246. Haeberli, W., Hoelzle, M., Paul, F., & Zemp, M. (2007). Integrated monitoring of mountain glaciers as key indicators of global climate change: the European Alps. Annals of Glaciology, 46, 150−160. Harper, J. T., Humphrey, N. F., Pfeffer, W. T., & Lazar, B. (2007). Two modes of accelerated glacier sliding related to water. Geophysical Research Letters, 34, L12503. doi:10.1029/ 2007GL030233 Hewitt, K. (2005). The Karakoram anomaly? Glacier expansion and the ‘elevation effect’, Karakoram Himalaya. Mountain Research and Development, 25(4), 332−340. Hubbard, B., & Glasser, N. F. (2005). Field Techniques in Glaciology and Glacial Geomorphology: Glacier Mass Balance and Motion (pp. 179−216). Chichester: Wiley. IPCC (2007a). Climate Change 2007: The Physical Science Basis. In S. Solomon, D. Qin, M. Manning, Z. Chen, M. Marquis, K. B. Averyt, M. Tignor, & H. L. Miller (Eds.), Contribution of Working Group I to the Fourth Assessment Report of the Intergovernmental Panel on Climate Change Cambridge: Cambridge University Press. IPCC (2007b). Climate Change 2007: Impacts, Adaption and Vulnerability. In M. L. Parry, O. F. Canziani, J. P. Palutikof, P. J. van der Linden, & C. E. Hanson (Eds.), Contribution of Working Group II to the Fourth Assessment Report of the Intergovernmental Panel on Climate Change Cambridge: Cambridge University Press. Joughin, I. (2002). Ice-sheet velocity mapping: a combined interferometric and speckletracking approach. Annals of Glaciology, 34, 195−201. Joughin, I., Abdalati, W., & Fahnestock, M. (2004). Large fluctuations in speed on Greenland's Jakobshavn Isbrae glacier. Nature, 432(7017), 608−610. Joughin, I., Kwok, R., & Fahnestock, M. (1996). Estimation of ice-sheet motion using satellite radar interferometry: method and error analysis with application to Humboldt Glacier, Greenland. Journal of Glaciology, 42(142), 564−575. Joughin, I., Slawek, T., Bindschadler, R., & Price, S. F. (2002). Changes in west Antarctic ice stream velocities: observation and analysis. Journal of Geophysical Research, 107 (B11), 2289. doi:10.1029/2001JB001029 Kääb, A. (2002). Monitoring high-mountain terrain deformation from repeated air- and spaceborne optical data: examples using digital aerial imagery and ASTER data. ISPRS Journal of Photogrammetry & Remote Sensing, 57, 39−52. Kääb, A. (2005). Combination of SRTM3 and repeat ASTER data for deriving alpine glacier flow velocities in the Bhutan Himalaya. Remote Sensing of Environment, 94, 463−474. Kargel, J. S., Abrams, M. J., Bishop, M. P., Bush, A., Hamilton, G., Jiskoot, H., et al. (2005). Multispectral imaging contributions to global land ice measurements from space. Remote Sensing of Environment, 99, 187−219. Karim, A., & Veizer, J. (2002). Water balance of the Indus River Basin and moisture source in the Karakoram and western Himalaya: Implications from hydrogen and oxygen isotopes in river water. Journal of Geophysical Research, 107(D18), 4362. doi:10.1029/2000JD000253 Kattelmann, R. (2003). Glacial lake outburst floods in the Nepal Himalaya: a manageable hazard? Natural Hazards, 28, 145−154. Khalsa, S. J. S., Dyurgerov, M. B., Khromova, T., Raup, B. H., & Barry, R. G. (2004). Spacebased mapping of glacier changes using ASTER and GIS tools. IEEE Transactions on Geoscience and Remote Sensing, 42(19), 2177−2183. Kulkarni, A. V., Rathore, B. P., Mahajan, S., & Mathur, P. (2005). Alarming retreat of Parbati glacier, Beas basin, Himachal Pradesh. Current Science, 88(11), 1844−1850. D. Scherler et al. / Remote Sensing of Environment 112 (2008) 3806–3819 Leprince, S., Barbot, S., Ayoub, F., & Avouac, J. P. (2007). Automatic and precise orthorectification, coregistration, and subpixel correlation of satellite images, application to ground deformation measurements. IEEE Transactions on Geoscience and Remote Sensing, 45(6), 1529−1558. Leprince, S., Berthier, E., Ayoub, F., Delacourt, C., & Avouac, J. P. (2008). Monitoring Earth surface dynamics with optical imagery. EOS, Transactions, 89(1), 1−2. Leprince, S., Musé, P., Avouac, J. P. (in press). In-Flight CCD Distortion Calibration for Pushbroom Satellites Based on Subpixel Correlation. IEEE Transactions on Geoscience and Remote Sensing, 46(9). doi:10.1109/TGRS.2008.918649 Liu, S., Shangguan, D., Ding, Y., Han, H., Xie, C., Zhang, Y., et al. (2006). Glacier changes during the past century in the Gangrigabo mountains, southeast Qinghai-Xizang (Tibetan) Plateau, China. Annals of Glaciology, 43, 187−193. Lucchita, B. K., & Ferguson, H. M. (1986). Antarctica: measuring glacier velocity from satellite images. Science, 234(4779), 1105−1108. Luckman, A., Quincey, D., & Bevan, S. (2007). The potential of satellite radar interferometry and feature tracking for monitoring flow rates of Himalayan glaciers. Remote Sensing of Environment, 111, 172−181. Luedeling, E., Siebert, S., & Buerkert, A. (2007). Filling the voids in the SRTM elevation model—A TIN-based delta surface approach. ISPRS Journal of Photogrammetry and Remote Sensing, 62(4), 283−294. Massonet, D., & Feigl, K. L. (1998). Radar interferometry and its application to changes in the Earth's surface. Reviews of Geophysics, 36(4), 441−500. Mayewski, P. A., & Jeschke, P. A. (1979). Himalayan and Trans-Himalayan glacier fluctuations since AD 1812. Arctic and Alpine Research, 11(3), 267−287. Michel, R., & Rignot, E. (1999). Flow of Glacier Moreno, Argentina, from repeat-pass shuttle imaging radar images: comparison of the phase correlation method with radar interferometry. Journal of Glaciology, 45(149), 93−100. Nye, J. F. (1959). A method of determining the strain-rate tensor at the surface of a glacier. Journal of Glaciology, 3(25), 409−419. Oerlemans, J. (1994). Quantifying global warming from the retreat of glaciers. Science, 264(5156), 243−245. Paul, F., Kääb, A., Maisch, M., Kellenberger, T., & Haeberli, W. (2002). The new remotesensing-derived Swiss glacier inventory: I. Methods. Annals of Glaciology, 34, 355−361. Rabatel, A., Dedieu, J. P., & Vincent, C. (2005). Using remote-sensing data to determine equilibrium-line altitude and mass-balance time series: validation on three French glaciers, 1994–2002. Journal of Glaciology, 51(175), 539−546. Reeh, N., Mohr, J. J., Madsen, S. N., Oerter, H., & Gundestrup, N. S. (2003). Threedimensional surface velocities of Storstrømmen glacier, Greenland, derived from radar interferometry and ice-sounding radar measurements. Journal of Glaciology, 49(165), 201−209. Rignot, E., & Kanagaratnam, P. (2006). Changes in the velocity structure of the Greenland Ice Sheet. Science, 311(5763), 986−990. 3819 Scambos, T. A., Dutkiewicz, M. J., Wilson, J. C., & Bindschadler, R. A. (1992). Application of image cross-correlation to the measurement of glacier velocity using satellite image data. Remote Sensing of Environment, 42, 177−186. Singh, P., Haritashya, U. K., & Kumar, N. (2007). Meteorological study for Gangotri Glacier and its comparison with other high altitude meteorological stations in central Himalayan region. Nordic Hydrology, 38(1), 59−77. Singh, P., Haritashya, U. K., Kumar, N., & Singh, Y. (2006). Hydrological characteristics of the Gangotri Glacier, central Himalayas, India. Journal of Hydrology, 327(1–2), 55−67. Stearns, L., & Hamilton, G. (2005). A new velocity map for Byrd Glacier, East Antarctica, from sequential ASTER satellite imagery. Annals of Glaciology, 41, 71−76. Strozzi, T., Luckman, A., Murray, T., Wegmüller, U., & Werner, C. L. (2002). Glacier motion estimation using SAR offset-tracking procedures. IEEE Transactions on Geoscience and Remote Sensing, 40(11), 2384−2391. Teshima, Y., & Iwasaki, A. (2008). Correction of attitude fluctuation terra spacecraft using ASTER/SWIR imagery with parallax observation. IEEE Transactions on Geoscience and Remote Sensing, 46(1), 222−227. Toutin, T. (2002a). Impact of terrain slope and aspect on radargrammetric DEM accuracy. ISPRS Journal of Photogrammetry & Remote Sensing, 57(3), 228−240. Toutin, T. (2002b). Three-dimensional topographic mapping with ASTER stereo data in rugged topography. IEEE Transactions on Geoscience and Remote Sensing, 40(10), 2241−2247. Toutin, T. (2004). Review article: geometric processing of remote sensing images: models, algorithms and methods. International Journal of Remote Sensing, 25(10), 1893−1924. Trouvé, E., Vasile, G., Gay, M., Bombrun, L., Grussenmeyer, P., Landes, T., et al. (2007). Combining airborne photographs and spaceborne SAR data to monitor temperate glaciers: potentials and limits. IEEE Transactions on Geoscience and Remote Sensing, 45(4), 905−924. Van Puymbroeck, N., Michel, R., Binet, R., Avouac, J. P., & Taboury, J. (2000). Measuring earthquakes from optical satellite images. Applied Optics, 39(20), 3486−3494. Wagnon, P., Linda, A., Arnaud, Y., Kumar, R., Sharma, P., Vincent, C., et al. (2007). Four years of mass balance on Chhota Shigri Glacier, Himachal Pradesh, India, a new benchmark glacier in the western Himalaya. Journal of Glaciology, 53(18), 603−611. Wessels, R. L., Kargel, J. S., & Kieffer, H. H. (2002). ASTER measurement of supraglacial lakes in the Mount Everest region of the Himalaya. Annals of Glaciology, 34, 399−408. Willis, I. C. (1995). Intra-annual variations in glacier motion: a review. Progress in Physical Geography, 19(1), 61−106. Winiger, M., Gumpert, M., & Yamout, H. (2005). Karakorum–Hindukush-western Himalaya: assessing high-altitude water resources. Hydrological Processes, 19, 2329−2338.