Pléiades Tri-Stereo Data for Glacier Investigations—Examples from the European Alps and the Khumbu Himal

In this study, we use Pléiades tri-stereo data to generate a digital elevation model (DEM) from the Pléiades images using a workflow employing semi-global matching (SGM). We examine the DEM accuracy in complex mountain glaciated terrain by comparing the new DEMs with an independent high-quality DEM based on airborne laser scanning (ALS) data for a study area in the Austrian Alps, and with ground control points for a study area in the Khumbu Himal of Nepal. The DEMs derived using the SGM algorithm compare well to the independent high-quality ALS DEM, and the workflow produces models of sufficient quality to resolve ground control points, which are based on Pléiades imagery that are of sufficient quality to perform high spatio-temporal resolution assessments of remote areas for which no field data is available. The relative accuracy is sufficient to investigate glacier surface elevation changes below one meter, and can therefore be applied over relatively short periods of time, such as those required for annual and seasonal assessments of change. The annual geodetic mass balance for the Alpine case derived from our DEM compares well to the glaciological mass balance, and multitemporal DEM analysis is used to resolve the seasonal changes of five glaciers in the Khumbu Himal, revealing that glaciological processes such as accumulation, ablation, and glacier movement mainly take place during the summer season, with the winter season being largely inactive in the year sampled.


Introduction
In recent years, the increasing availability and quality of satellite data has stimulated an expansion in the applications of satellite remote sensing to glaciological problems.Applications of satellite optical imagery in glaciology can be split into three main categories: (i) the derivation of digital elevation models (DEMs) from which volume changes and mass movements can be determined via DEM differencing, (ii) the identification of surface features based on optical properties, and (iii) the computation of glacier movement rates using image matching.
Glaciology has much to gain from satellite remote sensing, as glaciers are typically in remote and inaccessible locations for which the direct collection of data is logistically challenging.Of particular interest is the acquisition of glaciological information from poorly understood parts of the mountain cryosphere such as the mountains of high Asia where DEMs constructed from satellite imagery (e.g., ASTER, SRTM, SPOT5, Indian Remote Sensing Satellite (IRS-1C), or CORONA) are increasingly being used for glaciological studies.The most readily available sources of remote sensing-derived elevation data over the Himalayas is currently still the SRTM DEM version 3, but ASTER and Corona DEMs have been used in recent glaciological studies in the Nepal and Bhutan Himalaya, e.g., Bolch et al. [1] and Nuimura et al. [2].DEMs from the Indian Remote Sensing Satellite (IRS-1C) exist from a few areas such as the Baspa valley in the Himachal Pradesh district of India [3]; however, they are not in the public domain.Several studies focused on evaluating DEMs from stereo satellite data: Berthier et al. [4] compared SRTM elevations with SPOT5-derived elevations on non-glaciated terrain, while Bélart et al. [5] used WorldView and Pléiades data to investigate the winter mass balance of the Drangajökull ice cap in Iceland.The multi-annual surface lowering of glaciers obtained from DEM differencing has quantified the scale of the glacier mass losses in the Khumbu Himal [2, 6,7].Applications of surface classifications have included the detection and assessment of hazardous glacial lakes [8,9], the evaluation of ice avalanche potential in periglacial and glacial environments [10] or delimiting and mapping debris-covered glacier surfaces [3,[11][12][13].Developments in optical satellite remote sensing have led to a limited number of very-high resolution sensors, represented by IKONOS, QuickBird, OrbView-3, SPOT-5, Worldview 1-4, GeoEye, and the Pléiades twins.
In this study, we evaluate the use of very high-resolution (VHR) tri-stereo optical imagery (0.5-m resolution of panchromatic images) from the Pléiades-1 satellite constellation for glaciological applications.We apply photogrammetric techniques to Pléiades tri-stereo data and generate a 1-m resolution DEM of the study areas.Thereby, the stereoscopic triplet is in an enhanced acquisition mode compared to that of previous optical sensors (e.g., Quickbird, IKONOS).The tri-stereo images are composed of three nearly simultaneously acquired images, one backward (B) looking, one forward (F) looking, plus a third near-nadir (N) image [14].This configuration, in particular the use of the near-nadir image, allows a better retrieval of heights over terrains where the performance of classic photogrammetry with forward and backward-looking stereo pairs may be limited (e.g., high roughness, steep slopes, and shadows), and allows the calculation of dense photogrammetric point clouds by using three image pairs.The Pléiades mission is programmed to continue for several years [15], and the tri-stereo acquisition mode is set to become a common feature in future satellite missions.Recently, Pléiades data have been employed for DEM generation over mountainous terrain in several studies, estimating submeter to meter scale elevation changes that are associated with natural phenomena, such as earthquakes and glaciers [16].
Wagnon et al. [17] used SPOT5 and Pléiades DEM to calculate geodetic mass balance on Mera and Pokalde glacier (Nepal).Berthier et al. [16] tested the use of Pléiades imagery in glaciology on four different sites worldwide with promising results.Holzer et al. [18] used Pléiades DEM for glacier mass balance estimates at Muztag Ata (Pamir).Ruiz et al. [19,20] applied digital image matching to Pléiades data to generate accurate, high-resolution surface velocity maps, and used derived DEMs to calculate the glacier mass balances of the Monte Tronador glaciers in the Patagonian Andes.
An advantage of the Pléiades satellites to applications in complex, high-relief terrain is the ability to acquire tri-stereo data in a single overpass, which minimizes the problems caused by shadowing [21,22].In contrast to the previous studies, we employ redesigned pixel-matching algorithms and dense matching algorithms inherited from photogrammetry to generate DEMs from the Pleiades imagery [23].This paper presents an analysis of the quality and suitability of Pléiades tri-stereo data, which was processed using these photogrammetric techniques, for glaciological applications in mountainous terrain.First, the quality of Pléiades-derived DEM is evaluated against high accuracy airborne laser scanning data (ALS) in a well-studied region of the European Alps with relatively simple geometry and second, the same photogrammetric digital terrain model (DTM) generation workflow is applied to study glacier volume changes in a more complex terrain and the glacier geometries of the eastern Nepalese Himalaya.In addition, as has been demonstrated already by Ruiz et al. [19], image correlation is used to determine the glacier movement based on multitemporal orthoimages derived from the Pléiades data.

Study Areas
This study encompasses two different areas to explore the use of VHR remote sensing data to examine glacier behavior.A part of the Ötz valley in the Austrian Alps, in which substantial glaciological and environmental field data collected over the last six decades is available, as well as high-resolution remotely sensed data from airborne platforms, and the Khumbu Himal in the Himalayas of Nepal, for which field data is sparse and glacier behavior not well documented.
The study site in the Alps offers a wide array of available data, which is used to investigate the quality of our results.The high Himalaya provides more interesting challenges and unique opportunities to test and apply new area-wide remote sensing tools, such as the Pléiades twins, in order to improve the knowledge and understanding of the mountain cryosphere of the Himalaya.
The study area in the European Alps encompasses the Hochjochferner in the upper part of the Rofen valley, which is a side valley of the upper Ötz valley (see Figure 1).Hochjochferner is a relatively low-angled, predominantly clean-ice glacier that since the Little Ice Age has disintegrated into four separate glaciers, one of which was investigated in this study.In recent years, temperatures in the area have been above the long-term climatological mean [24], and Hochjochferner has been experiencing significant losses in area and mass during this time [25,26].Hochjochferner is investigated to compare the DEM calculated on the basis of Pléiades data (PLE-DEM) to a high-accuracy DEM based on airborne light detection and ranging (LiDAR) data (ALS-DEM), which is available for the area (see Section 3.1).

Study Areas
This study encompasses two different areas to explore the use of VHR remote sensing data to examine glacier behavior.A part of the Ötz valley in the Austrian Alps, in which substantial glaciological and environmental field data collected over the last six decades is available, as well as high-resolution remotely sensed data from airborne platforms, and the Khumbu Himal in the Himalayas of Nepal, for which field data is sparse and glacier behavior not well documented.
The study site in the Alps offers a wide array of available data, which is used to investigate the quality of our results.The high Himalaya provides more interesting challenges and unique opportunities to test and apply new area-wide remote sensing tools, such as the Pléiades twins, in order to improve the knowledge and understanding of the mountain cryosphere of the Himalaya.
The study area in the European Alps encompasses the Hochjochferner in the upper part of the Rofen valley, which is a side valley of the upper Ötz valley (see Figure 1).Hochjochferner is a relatively low-angled, predominantly clean-ice glacier that since the Little Ice Age has disintegrated into four separate glaciers, one of which was investigated in this study.In recent years, temperatures in the area have been above the long-term climatological mean [24], and Hochjochferner has been experiencing significant losses in area and mass during this time [25,26].Hochjochferner is investigated to compare the DEM calculated on the basis of Pléiades data (PLE-DEM) to a high-accuracy DEM based on airborne light detection and ranging (LiDAR) data (ALS-DEM), which is available for the area (see Section 3.1).In the Khumbu Himal of Nepal, a variety of glacier sizes, types, dynamics, topography, and debris cover is present (Figure 2).The wide altitude range and the variability in debris cover make the Himalayan glaciers particularly sensitive to climate forcing [27,28].Most importantly, field-based measurements of mass balance in the Himalaya are sparse, and Himalayan glaciers are conspicuously absent from global mass balance records [29].About 60% of the glaciated area investigated in this study is debris covered, and, in the absence of in situ measurements, remotely In the Khumbu Himal of Nepal, a variety of glacier sizes, types, dynamics, topography, and debris cover is present (Figure 2).The wide altitude range and the variability in debris cover make the Himalayan glaciers particularly sensitive to climate forcing [27,28].Most importantly, field-based measurements of mass balance in the Himalaya are sparse, and Himalayan glaciers are conspicuously absent from global mass balance records [29].About 60% of the glaciated area investigated in this study is debris covered, and, in the absence of in situ measurements, remotely sensed geodetic mass balances based on DEM differencing offers the best means of measuring current change in these glaciers.The role of the different processes contributing to the mass change of debris-covered Himalayan glaciers is still under investigation [30][31][32].In particular, avalanched snow is thought to contribute significantly to glacier accumulation, and rapid melting at exposed ice faces within the debris-covered zone have been shown to significantly influence the total glacier ablation [31], but both of these processes remain poorly quantified at the glacier scale.Here, we investigate the capabilities of the DEM calculated from Pléiades imagery to address those questions by observing the sub-annual and inter-annual surface changes of five glaciers: the Ngozumpa, Khumbu, Nuptse, Lhotse Nup, and Lhotse glaciers.They all flow southward, and are heavily debris-covered in the lower parts.While Nuptse, Lhotse Nup, and Lhotse are of similar dimensions, with similar lengths and terminate roughly at the same elevation, the Ngozumpa glacier, which originates on the summit regions of Cho Oyu with elevations of above 8000 m, and the Khumbu glacier are significantly larger and more complex, dendritic glaciers.

Pléiades Data and Processing
The Pléaides system consists of a constellation of two satellites for VHR panchromatic (PA) and multispectral (XS) optical observation of the Earth's surface.Pleiades 1A (PHR 1A) was launched in December 2011, while the second satellite, Pleiades 1B (PHR 1B) was launched in December 2012.Both satellites are operating on the same, sun-synchronous orbit with 98.2° inclination and an offset of 180° from each other, which together with their ability to tilt dramatically allows a minimum revisit time of any point on the globe within 24 h.The PA and XS images are acquired simultaneously at a ground sampling distance of 0.5 m and 2 m, respectively [32], whereas the multispectral dataset can be pan-sharpened to a resolution of 0.5 m [14,21].The Pléiades system is In particular, avalanched snow is thought to contribute significantly to glacier accumulation, and rapid melting at exposed ice faces within the debris-covered zone have been shown to significantly influence the total glacier ablation [31], but both of these processes remain poorly quantified at the glacier scale.Here, we investigate the capabilities of the DEM calculated from Pléiades imagery to address those questions by observing the sub-annual and inter-annual surface changes of five glaciers: the Ngozumpa, Khumbu, Nuptse, Lhotse Nup, and Lhotse glaciers.They all flow southward, and are heavily debris-covered in the lower parts.While Nuptse, Lhotse Nup, and Lhotse are of similar dimensions, with similar lengths and terminate roughly at the same elevation, the Ngozumpa glacier, which originates on the summit regions of Cho Oyu with elevations of above 8000 m, and the Khumbu glacier are significantly larger and more complex, dendritic glaciers.

Pléiades Data and Processing
The Pléaides system consists of a constellation of two satellites for VHR panchromatic (PA) and multispectral (XS) optical observation of the Earth's surface.Pleiades 1A (PHR 1A) was launched in December 2011, while the second satellite, Pleiades 1B (PHR 1B) was launched in December 2012.Both satellites are operating on the same, sun-synchronous orbit with 98.2 • inclination and an offset of 180 • from each other, which together with their ability to tilt dramatically allows a minimum revisit time of any point on the globe within 24 h.The PA and XS images are acquired simultaneously at a ground sampling distance of 0.5 m and 2 m, respectively [32], whereas the multispectral dataset can be pan-sharpened to a resolution of 0.5 m [14,21].The Pléiades system is the first of its kind that is capable of acquiring three or more nearly synchronous images of the same area with a stereo angle varying between ~6• and ~28 • , due to the high agility of the satellites [22].
For this study, Pléiades tri-stereo datasets were available from two acquisitions in the European Alps and three acquisitions in the Khumbu Himal (see Table 1), sampling a period of a year at each site.While the study area in the European Alps is of a size that could be covered with one Pléiades scene, the Khumbu Himal study area needed several images per acquisition.This poses particular challenges to the processing and analysis when different sectors of the study are covered by images from different dates (see Table 1 or Figures 3-5 for the exact acquisition dates), since, with acquisition times that are several weeks or even months apart, some surface changes are to be expected between the dates.Surface conditions, especially snow cover, can also be different within one acquisition period.This is clearly shown in the spring 2015 acquisition (see Figure 3), where the parts acquired on 8 April are still partly snow-covered, while there is no snow in the parts acquired on 20 May.
Remote Sens. 2017, 9, x FOR PEER REVIEW 5 of 26 the first of its kind that is capable of acquiring three or more nearly synchronous images of the same area with a stereo angle varying between ~6° and ~28°, due to the high agility of the satellites [22].For this study, Pléiades tri-stereo datasets were available from two acquisitions in the European Alps and three acquisitions in the Khumbu Himal (see Table 1), sampling a period of a year at each site.While the study area in the European Alps is of a size that could be covered with one Pléiades scene, the Khumbu Himal study area needed several images per acquisition.This poses particular challenges to the processing and analysis when different sectors of the study are covered by images from different dates (see Table 1 or Figures 3-5 for the exact acquisition dates), since, with acquisition times that are several weeks or even months apart, some surface changes are to be expected between the dates.Surface conditions, especially snow cover, can also be different within one acquisition period.This is clearly shown in the spring 2015 acquisition (see Figure 3), where the parts acquired on 8 April are still partly snow-covered, while there is no snow in the parts acquired on 20 May.For the spring 2015 and autumn 2015 acquisitions, three images were required to cover the Khumbu Himal study area, while four were needed for the spring 2016 acquisition (see Figure 5).Since the user does not have any influence on the tiling of the images in an acquisition, possible problems stemming from different distribution have to be kept in mind (see Section 4).For the spring 2015 and autumn 2015 acquisitions, three images were required to cover the Khumbu Himal study area, while four were needed for the spring 2016 acquisition (see Figure 5).Since the user does not have any influence on the tiling of the images in an acquisition, possible problems stemming from different distribution have to be kept in mind (see Section 4).DEM based on the Pléiades images were generated by first using semi-global matching (SGM) [33] within the IMAGINE photogrammetry suite of ERDAS IMAGINE.This results in a dense photogrammetric point cloud.Since the rational polynomial coefficient (RPC) files are delivered with the original scene, orientation and georectification for the Pléiades imagery can be achieved through established and well-documented processing e.g., Haala and Tiris [34,35].The RPC file that is provided with the Pléiades data is used to define the initial functions for transforming the sensor geometry to image geometry.With those transformation functions, the individual geometries of   DEM based on the Pléiades images were generated by first using semi-global matching (SGM) [33] within the IMAGINE photogrammetry suite of ERDAS IMAGINE.This results in a dense photogrammetric point cloud.Since the rational polynomial coefficient (RPC) files are delivered with the original scene, orientation and georectification for the Pléiades imagery can be achieved through established and well-documented processing e.g., Haala and Tiris [34,35].The RPC file that is provided with the Pléiades data is used to define the initial functions for transforming the sensor geometry to image geometry.With those transformation functions, the individual geometries of DEM based on the Pléiades images were generated by first using semi-global matching (SGM) [33] within the IMAGINE photogrammetry suite of ERDAS IMAGINE.This results in a dense photogrammetric point cloud.Since the rational polynomial coefficient (RPC) files are delivered with the original scene, orientation and georectification for the Pléiades imagery can be achieved through established and well-documented processing e.g., Haala and Tiris [34,35].The RPC file that is provided with the Pléiades data is used to define the initial functions for transforming the sensor geometry to image geometry.With those transformation functions, the individual geometries of each image in the triplet are then orientated relative to each other.To obtain the most accurate exterior orientation possible, the initial RPC functions were refined using tie points (TP).We used automatically identified and manually controlled TPs over the PA triplets, obtaining an overall root mean square error of less than 0.1 pixels (<0.05 m) shift between images at each TP within all of the acquisitions.For the data in the Hochjochferner area, 12 ground control points (GCPs) collected from a high-resolution orthoimage and an ALS DEM, which was available for the state of Tyrol from 2013 [36], were used in addition to the automatically collected TPs.GCPs for the Khumbu Himal study area were not available in a suitable spatial distribution, so the processing was solely based on the RPC files and TPs.
The SGM algorithm was used to calculate point clouds (x, y, z) from each possible image pair (F-B, F-N, and N-B), which resulted in three point clouds for each Pléiades triplet.For a detailed description of the generation of the photogrammetric point clouds, see Galos et al. [37].The calculated point clouds for the N-B, F-N, and F-B combinations were then filtered for outliers, which are mainly found in very steep and shaded areas, using a local topographic three-dimensional (3D) filter.By merging the three point clouds, we obtained a dataset with a higher point density and less gaps, since the three image pairs offer higher chances of successful cross-correlation, especially in steep or shaded areas.The merged point cloud was then gridded to a spatial resolution of 1 × 1 m pixel size using the average elevation of all of the points within one raster cell as the elevation value for the cell.
Gaps were still present in very steep areas as well as in areas with low contrast because of clouds, fresh snow, or liquid water.Gaps only occurred in about 2% of the study area, but are irregularly distributed.The main problem with the calculation of the photogrammetric point clouds, and respectively the DEM, occurred in areas with low contrast due to relatively fresh snow.Those are mainly located in the highest parts of the Ngozumpa glacier (see Figures 3-5).The DEM shows a lot of noise in those areas, which results in an unrealistic representation of the terrain surface (see Figure 6).To account for the flat-water surfaces of proglacial lakes, we mapped unfrozen water bodies from the optical images and set those areas as flat surfaces in the DEM by using the mean value of the surrounding cells.
Remote Sens. 2017, 9, x FOR PEER REVIEW 7 of 26 each image in the triplet are then orientated relative to each other.To obtain the most accurate exterior orientation possible, the initial RPC functions were refined using tie points (TP).We used automatically identified and manually controlled TPs over the PA triplets, obtaining an overall root mean square error of less than 0.1 pixels (<0.05 m) shift between images at each TP within all of the acquisitions.For the data in the Hochjochferner area, 12 ground control points (GCPs) collected from a high-resolution orthoimage and an ALS DEM, which was available for the state of Tyrol from 2013 [36], were used in addition to the automatically collected TPs.GCPs for the Khumbu Himal study area were not available in a suitable spatial distribution, so the processing was solely based on the RPC files and TPs.
The SGM algorithm was used to calculate point clouds (x, y, z) from each possible image pair (F-B, F-N, and N-B), which resulted in three point clouds for each Pléiades triplet.For a detailed description of the generation of the photogrammetric point clouds, see Galos et al. [37].The calculated point clouds for the N-B, F-N, and F-B combinations were then filtered for outliers, which are mainly found in very steep and shaded areas, using a local topographic three-dimensional (3D) filter.By merging the three point clouds, we obtained a dataset with a higher point density and less gaps, since the three image pairs offer higher chances of successful cross-correlation, especially in steep or shaded areas.The merged point cloud was then gridded to a spatial resolution of 1 × 1 m pixel size using the average elevation of all of the points within one raster cell as the elevation value for the cell.
Gaps were still present in very steep areas as well as in areas with low contrast because of clouds, fresh snow, or liquid water.Gaps only occurred in about 2% of the study area, but are irregularly distributed.The main problem with the calculation of the photogrammetric point clouds, and respectively the DEM, occurred in areas with low contrast due to relatively fresh snow.Those are mainly located in the highest parts of the Ngozumpa glacier (see Figures 3-5).The DEM shows a lot of noise in those areas, which results in an unrealistic representation of the terrain surface (see Figure 6).To account for the flat-water surfaces of proglacial lakes, we mapped unfrozen water bodies from the optical images and set those areas as flat surfaces in the DEM by using the mean value of the surrounding cells.

ALS Data and Processing
Airborne laser scanning data was available for the study area in the European Alps from 2013.This was used to investigate the accuracy of the calculated DEM based on Pléiades imagery.The ALS data were acquired on 5 September 2013, with an average point density of more than 2.5 points per square meter.For more information on the ALS data, the reader is referred to Rieg et al. [38], where data from the same acquisition was used and is described in more detail.
The generation of the DEM based on the available ALS data was straightforward.Data was received as a filtered and georeferenced point cloud.Those points were imported into a Laserdata Information System Database [39], and a DEM with a cell size of 1 m was created by using the average elevation value of the points within a potential raster cell as the cell value.Due to the high point density of the available ALS data, the resulting DEM had gaps in less than 1% of the cells, which were filled using the mean value of the 20 closest cells containing data as the new cell value.This resulted in an accurate representation of the terrain and its features.The accuracy of this DEM was evaluated on the basis of a comparison to differential Global Navigation Satellite System (dGNSS) measurements in the study area and to a flat and even reference surface close to the study area [40,41].

Co-Registration
The final DEMs for each study site were co-registered using the method described by Berthier [42], which uses the relationship of surface differences in stable terrain with terrain slope and aspect to correct the horizontal shifts between different DEMs of the same area.The ALS-DEM from 2013 was used as the master DEM for the study area in the European Alps, and the most recent DEM (spring 2016) as the master DEM for the Khumbu Himal, as this coincided with the acquisition of seven GCPs on 12-13 April 2016, in the southwestern part of the Khumbu Himal study area around the Ngozumpa glacier terminus.The remaining "slave" DEMs were all tested against the respective master DEM, which means the differential DEM (dDEM) was related to terrain slope and aspect [42].The shift in the x, y, and z-direction of the slave DEM in relation to the master DEM (Table 2) was corrected as suggested by Berthier [42].

Analysis of Glacier Processes
To analyze glacier surface changes, dDEMs were calculated by subtracting the older DEM from the newer DEM of the respective period.Since gaps pose a problem for the analysis, they had to be filled.This was done in the dDEMs to avoid creating very smooth areas in the DEMs, which would lead to erroneous difference values in rugged mountain terrain.The gaps in the dDEMs were filled by using the mean value of the next 20 cells containing data.
Firstly, the outlines of Hochjochferner and the Khumbu Himal region were determined.On Hochjochferner, this was done by manually adjusting the data from the Austrian glacier inventory [41] based on the 2014 Pléiades orthoimages, while for the Khumbu Himal region, outlines for the five study glaciers (Figure 2) were manually mapped on the basis of the orthoimages and DEM from spring 2016.
Multitemporal orthophotos were generated from the available Pléiades scenes.Those can then be used to determine the glacier movement by using image correlation software.In our case, we used IMCORR [43] to obtain surface displacement rates.The IMCORR software determines the surface displacement rates by using a reference block of the raster cell values of the first image and correlating the values within a determined search block.The output consists of surface displacement vectors.

Comparison of Pléiades DEM with ALS-DEM in the European Alps
To investigate the absolute accuracy of the DEM calculated from Pléiades tri-stereo data in high mountain areas, the data from the Eastern Alps was compared with the ALS-DEM, which was taken to be 'correct'.Co-registration of the two Pleaides DEMs to the ALS-DEM revealed significant horizontal offsets and more minor vertical offsets (Table 2).Results can be substantially improved (Figure 7) by shifting the respective DEM according to these offsets, as suggested by Berthier [42].
Remote Sens. 2017, 9, x FOR PEER REVIEW 9 of 26 five study glaciers (Figure 2) were manually mapped on the basis of the orthoimages and DEM from spring 2016.Multitemporal orthophotos were generated from the available Pléiades scenes.Those can then be used to determine the glacier movement by using image correlation software.In our case, we used IMCORR [43] to obtain surface displacement rates.The IMCORR software determines the surface displacement rates by using a reference block of the raster cell values of the first image and correlating the values within a determined search block.The output consists of surface displacement vectors.

Comparison of Pléiades DEM with ALS-DEM in the European Alps
To investigate the absolute accuracy of the DEM calculated from Pléiades tri-stereo data in high mountain areas, the data from the Eastern Alps was compared with the ALS-DEM, which was taken to be 'correct'.Co-registration of the two Pleaides DEMs to the ALS-DEM revealed significant horizontal offsets and more minor vertical offsets (Table 2).Results can be substantially improved (Figure 7) by shifting the respective DEM according to these offsets, as suggested by Berthier [42].To assess the errors remaining after co-registration, the surface difference between each Pleiades-DEM and the ALS-DEM in the areas that were identified as likely to be stable was investigated.Although stable areas are scarce within the mountain environment, three areas with a total size of ca.15,000 m 2 were defined, with a mean absolute surface difference to the ALS-DEM of 0.78 m and 0.60 m, for the 2014 and 2015 Pléiades DEM, respectively.These submeter residual differences demonstrate that through using the processing methods applied here, it is possible generate DEMs from the tri-stereo Pléiades images that offer a sufficient accuracy for many applications, even in mountainous areas.Since the errors remaining after the co-registration are non-linear, they are not easily corrected.Especially in the northwestern part of the study area, a lower number of TP due to the snow-covered surface in all of the acquisitions can also lead to a decreased fit between the different DEM.3.2 relative orientation of Pléiades DEM in the Khumbu Himal study area.
To cope with the horizontal inaccuracies in the Khumbu Himal data, the DEM of spring and autumn 2015 in the Khumbu Himal were co-registered to the 2016 data (see Section 2.2.3 and Table 2).To better understand the remaining vertical error after the co-registration, 12 stable areas with inclination <40° and a total size of 6.5 km 2 were selected throughout the study area.The mean absolute difference of the two 2015 DEM to the 2016 DEM was calculated, and shows an average of 0.60 m for spring 2015 and 0.40 m for autumn 2015, and the standard deviation of the surface elevation difference in the stable areas for all three differences is shown in Table 3.This suggests that while the co-registration to the spring 2016 DEM does not necessarily increase the absolute accuracy, as we do not know which of the DEMs is most accurate, the relative accuracy within our dataset is To assess the errors remaining after co-registration, the surface difference between each Pleiades-DEM and the ALS-DEM in the areas that were identified as likely to be stable was investigated.Although stable areas are scarce within the mountain environment, three areas with a total size of ca.15,000 m 2 were defined, with a mean absolute surface difference to the ALS-DEM of 0.78 m and 0.60 m, for the 2014 and 2015 Pléiades DEM, respectively.These submeter residual differences demonstrate that through using the processing methods applied here, it is possible generate DEMs from the tri-stereo Pléiades images that offer a sufficient accuracy for many applications, even in mountainous areas.Since the errors remaining after the co-registration are non-linear, they are not easily corrected.Especially in the northwestern part of the study area, a lower number of TP due to the snow-covered surface in all of the acquisitions can also lead to a decreased fit between the different DEM.3.2 relative orientation of Pléiades DEM in the Khumbu Himal study area.
To cope with the horizontal inaccuracies in the Khumbu Himal data, the DEM of spring and autumn 2015 in the Khumbu Himal were co-registered to the 2016 data (see Section 2.2.3 and Table 2).To better understand the remaining vertical error after the co-registration, 12 stable areas with inclination <40 • and a total size of 6.5 km 2 were selected throughout the study area.The mean absolute difference of the two 2015 DEM to the 2016 DEM was calculated, and shows an average of 0.60 m for spring 2015 and 0.40 m for autumn 2015, and the standard deviation of the surface elevation difference in the stable areas for all three differences is shown in Table 3.This suggests that while the co-registration to the spring 2016 DEM does not necessarily increase the absolute accuracy, as we do not know which of the DEMs is most accurate, the relative accuracy within our dataset is sufficient to permit DEM differencing analysis for processes resulting in surface elevation changes >1 m after performing co-registration.The surface difference in off-glacier terrain has been investigated to give an estimation of the error within the surface elevation differences.The surface differences along two profiles close to the glacier tongues of the Ngozumpa and Khumbu glacier show the errors in the dDEM (Figures 8 and 9), while their locations are shown in Figure 10.Both profiles show differences of up to ±1 m in stable areas, with very few outliers of up to ±7 m in profile 2.
Remote Sens. 2017, 9, x FOR PEER REVIEW 10 of 26 sufficient to permit DEM differencing analysis for processes resulting in surface elevation changes >1 m after performing co-registration.The surface difference in off-glacier terrain has been investigated to give an estimation of the error within the surface elevation differences.The surface differences along two profiles close to the glacier tongues of the Ngozumpa and Khumbu glacier show the errors in the dDEM (Figures 8 and  9), while their locations are shown in Figure 10.Both profiles show differences of up to ±1 m in stable areas, with very few outliers of up to ±7 m in profile 2.   sufficient to permit DEM differencing analysis for processes resulting in surface elevation changes >1 m after performing co-registration.The surface difference in off-glacier terrain has been investigated to give an estimation of the error within the surface elevation differences.The surface differences along two profiles close to the glacier tongues of the Ngozumpa and Khumbu glacier show the errors in the dDEM (Figures 8 and  9), while their locations are shown in Figure 10.Both profiles show differences of up to ±1 m in stable areas, with very few outliers of up to ±7 m in profile 2.

Accuracy of Khumbu Himal Pléiades DEM based on RPC
At the moment, no other high-accuracy DEM of the Khumbu Himal region is publicly available to assess both the accuracy of the heights and their horizontal position in the derived Pléiades DEM.We therefore used differential GNSS solutions for the x, y, and z positions of seven GCPs that were measured in April 2016.The corners of four roofs in the villages of Gokyo and Dragnag, as well as at a ruin south of Gokyo, and two massive boulders were measured with differential GNSS using a Trimble Geo7X with Zephyr antenna and a local base station.Point positions were collected at 1-s intervals over a period of at least 10 min and post-processed differentially with the base station to give a high precision of the location.

Accuracy of Khumbu Himal Pléiades DEM based on RPC
At the moment, no other high-accuracy DEM of the Khumbu Himal region is publicly available to assess both the accuracy of the heights and their horizontal position in the derived Pléiades DEM.We therefore used differential GNSS solutions for the x, y, and z positions of seven GCPs that were measured in April 2016.The corners of four roofs in the villages of Gokyo and Dragnag, as well as at a ruin south of Gokyo, and two massive boulders were measured with differential GNSS using a Trimble Geo7X with Zephyr antenna and a local base station.Point positions were collected at 1-s intervals over a period of at least 10 min and post-processed differentially with the base station to give a high precision of the location.
Due to the uneven distribution of GCPs within the study area, these points were not incorporated into the DEM generation, but are instead used as a measure to check the accuracy of the spring 2016 DEM used as the master DEM in the Khumbu Himal study area, which was based on imagery acquired three weeks earlier than the GCPs.
The horizontal accuracy was estimated from offsets between the x and y coordinates of the GCPs and their clearly visible positions in a pan-sharpened image that was formed using the near-nadir PA/XS images and orthorectified using the DEM processed solely with the use of RPC information (Table 4).It should be noted that the distances were measured from the corrected dGNSS position to the center of the cell defined as the respective corner.Therefore, the offset is not necessarily the true distance, but it still gives an overview of the accuracy of the DEM at its calculated resolution.Due to the uneven distribution of GCPs within the study area, these points were not incorporated into the DEM generation, but are instead used as a measure to check the accuracy of the spring 2016 DEM used as the master DEM in the Khumbu Himal study area, which was based on imagery acquired three weeks earlier than the GCPs.
The horizontal accuracy was estimated from offsets between the x and y coordinates of the GCPs and their clearly visible positions in a pan-sharpened image that was formed using the near-nadir PA/XS images and orthorectified using the DEM processed solely with the use of RPC information (Table 4).It should be noted that the distances were measured from the corrected dGNSS position to the center of the cell defined as the respective corner.Therefore, the offset is not necessarily the true distance, but it still gives an overview of the accuracy of the DEM at its calculated resolution.The combined horizontal offset reaches about 2.5 m in some cases, which is in line with the results from previous studies [44] and significantly lower than the 8-10 m offsets reported in other studies [41], which can be considered a satisfying result from the RPC, TP, and co-registration workflow.The comparison also shows there to be no systematic offset between the DEM and the GCPs, as even the values for neighbouring GCP locations (e.g., Gokyo and Gokyo hospital or Dragnag 1 and Dragnag 2) show offsets in different directions.Therefore, as was also noted by Perko et al. [33], an improvement of the geolocation through a pixel shift to compensate a systematic locational bias is not possible in this case.The vertical difference between the DEM and the dGNSS measurements are larger than the horizontal ones (see Table 3).However, since the relative orientation of the DEMs is better than the deviation of the DEM in relation to the dGNSS measurements, this is of little consequence for the analysis of surface elevation change.Nevertheless, this rather low vertical accuracy of DEM calculated without GCP from Pléiades data has to be kept in mind if it is used for analyses in conjunction with other data sources, such as TanDEM-X, SRTM, and ASTERGDEM.

Annual Surface Elevation Change at Hochjochferner
Surface elevation change at Hochjochferner is negative in all elevation bands over the mass balance year studied, though only slightly so in the highest elevation bands, and reaching a maximum at the glacier terminus (see Figures 11 and 12).Therefore, the glacier is lowering over its whole area, indicating a glacier that is strongly out of balance with the prevailing climate of this year.Only very small areas in the southwestern part of the upper glacier body show slight increases in surface elevation, which is interpreted as a localized accumulation of snow through the avalanches in this area.This agrees with the glaciological measurements of climatological surface mass balance on this glacier over the same period, which identified only small, localized areas of mass gain, and a specific glacier mass balance of −2030 mm water equivalent [45].

Surface Elevation Changes in the Khumbu Himal
Since the large Khumbu and Ngozumpa glaciers are located in two parts of the Pléiades acquisitions (see Figures 3-5), it was necessary to create mosaics based on the individual Pléiades scenes of the Khumbu Himal study area.While this is a straightforward process, the fit between the individual parts within the final mosaic was impacted by the DEM being calculated without the use of GCPs.The difference in z-values for the overlapping areas ranged from 0.3 m for the southwestern and northwestern DEM for the spring 2016 acquisition to 1.2 m for the southwestern and northwestern DEM for the spring 2015 acquisition.Due to fringe effects, which make the fringes of the calculated DEM not usable, the overlapping areas of the different images were too small to properly co-register the parts of the different DEM of a single acquisition.Therefore, we did not use data based on mosaicked DEM, but mosaicked the calculated dDEM.This has the advantage that, by solely using the differences, which are calculated based on already well-co-registered DEM, the non-perfect fit of the different DEM of one acquisition to each other is less of an issue.For the overlapping areas of the respective acquisition, the mean value of the two dDEM covering the area was used in the final dDEM for the respective glacier.
Despite the efforts regarding the correction of distortion and misalignments of the DEM and dDEM, data from different images of one acquisition does not always fit perfectly in some areas of the Khumbu Himal study area.This is visible in the dDEM with surface differences of more than 0.5 m in stable areas, or steps in the surface elevation changes, e.g., on the tongue of Ngozumpa glacier for the winter surface elevation difference (see Figure 13).

Surface Elevation Changes in the Khumbu Himal
Since the large Khumbu and Ngozumpa glaciers are located in two parts of the Pléiades acquisitions (see Figures 3-5), it was necessary to create mosaics based on the individual Pléiades scenes of the Khumbu Himal study area.While this is a straightforward process, the fit between the individual parts within the final mosaic was impacted by the DEM being calculated without the use of GCPs.The difference in z-values for the overlapping areas ranged from 0.3 m for the southwestern and northwestern DEM for the spring 2016 acquisition to 1.2 m for the southwestern and northwestern DEM for the spring 2015 acquisition.Due to fringe effects, which make the fringes of the calculated DEM not usable, the overlapping areas of the different images were too small to properly co-register the parts of the different DEM of a single acquisition.Therefore, we did not use data based on mosaicked DEM, but mosaicked the calculated dDEM.This has the advantage that, by solely using the differences, which are calculated based on already well-co-registered DEM, the non-perfect fit of the different DEM of one acquisition to each other is less of an issue.For the overlapping areas of the respective acquisition, the mean value of the two dDEM covering the area was used in the final dDEM for the respective glacier.
Despite the efforts regarding the correction of distortion and misalignments of the DEM and dDEM, data from different images of one acquisition does not always fit perfectly in some areas of the Khumbu Himal study area.This is visible in the dDEM with surface differences of more than 0.5 m in stable areas, or steps in the surface elevation changes, e.g., on the tongue of Ngozumpa glacier for the winter surface elevation difference (see Figure 13).Some noteworthy glaciological features are revealed by the DEM differencing in the Khumbu Himal study area.Since the glaciers in the Khumbu Himal are of the summer accumulation type [46,47], surface elevation change rates are very low in winter.Very little change on the glaciers seems to occur between October and April (see Figures 13-16), although there seems to be some winter accumulation in the highest parts of the Lhotse Nup glacier, which fits to the findings of [48].As surface elevation change mainly occurs in summer, the altitudinal distribution of the surface elevation change over the whole year is similar to the one just for the summer period (see Figures 13-16).The lower areas of all of the glaciers show very little average surface elevation change, although there seems to be a slight positive change for the Lhotse Nup and Khumbu glaciers, and a slight decrease in the average surface elevation for the Lhotse and Ngozumpa glaciers.Decreases in surface elevation are mainly located in the lower elevation bands, which also have larger areas, while surface increase occurs in the higher elevation bands with rather small areas.The total volume change of the four glaciers is shown in Table 5.Some noteworthy glaciological features are revealed by the DEM differencing in the Khumbu Himal study area.Since the glaciers in the Khumbu Himal are of the summer accumulation type [46,47], surface elevation change rates are very low in winter.Very little change on the glaciers seems to occur between October and April (see Figures 13-16), although there seems to be some winter accumulation in the highest parts of the Lhotse Nup glacier, which fits to the findings of [48].As surface elevation change mainly occurs in summer, the altitudinal distribution of the surface elevation change over the whole year is similar to the one just for the summer period (see Figures 13-16).The lower areas of all of the glaciers show very little average surface elevation change, although there seems to be a slight positive change for the Lhotse Nup and Khumbu glaciers, and a slight decrease in the average surface elevation for the Lhotse and Ngozumpa glaciers.Decreases in surface elevation are mainly located in the lower elevation bands, which also have larger areas, while surface increase occurs in the higher elevation bands with rather small areas.The total volume change of the four glaciers is shown in Table 5.Secondly, surface elevation change rates on the heavily debris-covered glacier tongues are relatively low, with the highest changes at exposed ice cliffs.Nearly no winter season surface displacement is visible on the glacier tongues (see Figures 13 and 16).
The results for the highest parts of the Ngozumpa glacier, which are located above 7000 m a.s.l., show high changes in surface elevation, especially between the two spring acquisitions and the spring 2015 and autumn 2015 acquisitions (see Figures 17 and 18     The results for the highest parts of the Ngozumpa glacier, which are located above 7000 m.a.s.l., show high changes in surface elevation, especially between the two spring acquisitions and the spring 2015 and autumn 2015 acquisitions (see Figures 17 and 18).The values, which show surface elevation changes of more than 15 m over large areas, are unrealistic and most likely caused by errors in the DEM.Those areas, while shown in Figures 18 and 19 for reference, are therefore not included in the following analysis.

Image Correlation for Deriving Glacier Surface Velocity
Since large parts of the investigated glaciers were snow-covered in the spring 2015 acquisition, which prohibits a sufficient correlation of the image blocks, only the panchromatic orthoimages from autumn 2015 and spring 2016 with a resolution of 0.5 m were used as input data for IMCORR (see Section 2.2.4).Even for this image pair, the glaciers in the eastern part of the study area showed a significant snow cover in the higher parts, which made the image matching in these areas difficult, which is why we only present results for the Ngozumpa glacier (see Figure 17).Movement rates

Image Correlation for Deriving Glacier Surface Velocity
Since large parts of the investigated glaciers were snow-covered in the spring 2015 acquisition, which prohibits a sufficient correlation of the image blocks, only the panchromatic orthoimages from autumn 2015 and spring 2016 with a resolution of 0.5 m were used as input data for IMCORR (see Section 2.2.4).Even for this image pair, the glaciers in the eastern part of the study area showed a significant snow cover in the higher parts, which made the image matching in these areas difficult, which is why we only present results for the Ngozumpa glacier (see Figure 17).Movement rates increase in the higher parts of the Ngozumpa glacier, with surface displacement exceeding 10 m over the winter investigation period in some areas.The results also clearly show which tributary areas are still dynamically connected to the main glacier and the areas that conjoin the glacier, but do not seem to be part of its dynamics anymore [49].This is especially evident in the case of some of the lower branches, which do not show any movement; this could connect them with the main glacier tongue.

Problems of DEM based on Pléiades Data in High Mountain Areas
A general problem of optical remote sensing data is the possibility of cloud cover during the time of acquisition.It is possible to specify a cloud cover of less than 5% of the area of interest for the Pléaides acquisition, but of course, even a low percentage of cloud cover can obscure important parts

Problems of DEM based on Pléiades Data in High Mountain Areas
A general problem of optical remote sensing data is the possibility of cloud cover during the time of acquisition.It is possible to specify a cloud cover of less than 5% of the area of interest for the Pléaides acquisition, but of course, even a low percentage of cloud cover can obscure important parts in the acquired images, e.g., the location of ground control points or objects that are to be studied.A low accepted percentage of cloud cover for a planned acquisition will also lower the chance of an actual acquisition taking place, since the operating company of the satellites will only acquire images if the feasibility of meeting the parameters of the order (including cloud cover) is high.While this is less problematic in areas with rather predictable weather conditions, such as the Nepalese Himalayas with their seasonal climate, it is harder in the Alps or in tropical areas with frequent cloud cover throughout the year.Nevertheless, in the Khumbu Himal study site, meeting the cloud cover criteria resulted in seasonal acquisitions that spanned several dates, which brings data processing and mosaicking challenges.
Fresh snow can also be problematic for the calculation of DEM based on optical stereo data in several ways.The low contrast of the white snow surfaces leads to fewer matches being identified within the matching algorithm and therefore more gaps, outliers, and higher noise levels in the data over bright snow-covered areas.Of course, with snow covering the ground, the snow surface will be represented in the DEM, and this will be located above the real terrain surface as the snow fills terrain structures such as holes or furrows.This can lead to significant differences in rugged mountain areas, where large amounts of snow can accumulate in some areas due to topography and wind drift.This makes it crucial that data is acquired during minimum snow-cover conditions if the goal of the study is the accurate detection of glacier surface elevation changes and the correct calculation of volume changes of glaciers.
Unfrozen water bodies are also problematic, since their surface changes between the taking of the stereo images, which is problematic for the matching algorithm.This leads to either missing or an erroneous elevation values for water bodies.Since there are many ponds on the heavily debris-covered tongues of the glaciers in the Khumbu Himal, this has to be considered.

Accuracy of DEM based on Pléiades Tri-Stereo Data With and Without the Use of GCPs
The DEM for the study area in the Alps, which were calculated using GCPs, show horizontal offsets of similar magnitude to the ones for the Khumbu Himal, which were calculated without GCP (see Tables 2 and 3), suggesting that horizontal accuracy for DEMs calculated from Pléiades tri-stereo images does not improve significantly if GCPs are included in the DEM generation.However, using GCPs appears to substantially improve the vertical accuracies.This agrees with the findings of Salerno et al. [48], who reported that GCPs mainly improved the vertical accuracy.
As another measure of the accuracy of the calculated results, and correspondingly the applicability of the methods, the different periods of available surface difference data for the Khumbu Himal study area are compared.That the combined surface elevation changes of the winter and summer seasons compare well to the annual surface elevation changes between spring 2015 and spring 2016 provide evidence that the dDEM data is of sufficient quality to resolve sub-annual processes on Himalayan glaciers (see Figure 19).
This comparison shows that the results are similar for the three rather similar glaciers in the Khumbu Himal (Lhotse Nup, Nuptse, and Lhotse glaciers, see Figure 19).Only in the lower parts of Nuptse glacier are some rather small differences are visible in several elevation bands.The larger and more complex Ngozumpa and Khumbu are behaving differently with slightly higher differences between the two curves.
The distribution over the elevation bands as well as the magnitude of the surface elevation change is nearly identical in most cases.We therefore conclude that differential DEMs from Pléiades tri-stereo data provide meaningful results, which can be seen as a proof of concept for the applicability of Pléiades data for surface elevation change calculations in areas without field data to constrain the DEM accuracy.The relative accuracy of the DEM is higher than expected, even though the vertical actual accuracy, e.g., in comparison to GCP measured in the area, is lower than reported by previous studies.If this, and that problems which can arise in very steep and snow-covered terrain are kept in mind, Pléiades DEM are well-suited for a variety of applications, even in areas where no GCPs in sufficient number, accuracy, or spatial distribution are available.

Glaciological Interpretation
The glacier-wide surface lowering observed at Hochjochferner is in good accordance with the results from other glaciers in the area [25,26] and the glaciers of the adjacent Vinschgau catchment [38].Surface lowering near the terminus reaches the order of 5 m per year.
The results for the glaciers in the Khumbu Himal study area differ significantly from those in the Alpine case.The typical surface lowering of the debris-covered surface is much less than in the European Alps, at least in relation to their large areas.An increase in surface elevation is concentrated in the highest parts of the glaciers and restricted to the areas identified as avalanche cones in the optical imagery (see Figure 20).In the studied year, this avalanche accumulation mainly occurs only during the summer period, while areas showing a significant increase in surface elevation during winter are rare.

Glaciological Interpretation
The glacier-wide surface lowering observed at Hochjochferner is in good accordance with the results from other glaciers in the area [25,26] and the glaciers of the adjacent Vinschgau catchment [38].Surface lowering near the terminus reaches the order of 5 m per year.
The results for the glaciers in the Khumbu Himal study area differ significantly from those in the Alpine case.The typical surface lowering of the debris-covered surface is much less than in the European Alps, at least in relation to their large areas.An increase in surface elevation is concentrated in the highest parts of the glaciers and restricted to the areas identified as avalanche cones in the optical imagery (see Figure 20).In the studied year, this avalanche accumulation mainly occurs only during the summer period, while areas showing a significant increase in surface elevation during winter are rare.This is despite occasional winter storms being reported for the region.A decrease in surface elevation is mainly a summer phenomenon, and it is concentrated in the middle part of the glacier's long profiles, contributing to the development of the low-angled glacier termini and surface water storage that is typical of debris-covered glaciers in the region [50].A part of the changes seems to occur at exposed ice cliffs, while the debris-covered part of the glaciers show rather little change, especially in the lower parts of the Lhotse, Nuptse, and Lhotse Nup glaciers, while there is some negative surface elevation change on the tongues of the Ngozumpa and Khumbu glaciers.This is despite occasional winter storms being reported for the region.A decrease in surface elevation is mainly a summer phenomenon, and it is concentrated in the middle part of the glacier's long profiles, contributing to the development of the low-angled glacier termini and surface water storage that is typical of debris-covered glaciers in the region [50].A part of the changes seems to occur at exposed ice cliffs, while the debris-covered part of the glaciers show rather little change, especially in the lower parts of the Lhotse, Nuptse, and Lhotse Nup glaciers, while there is some negative surface elevation change on the tongues of the Ngozumpa and Khumbu glaciers.
Of course, one year of data for only a few glaciers cannot be used to draw general conclusions about the behavior of glaciers in the area, but the results indicate that net accumulation and ablation on the glaciers in the Khumbu Himal mainly occur during the summer.Accumulation seems to be limited to mass input by avalanches or at least strongly depends on it, and surface change over the large debris-covered parts of the glaciers is minimal.

Conclusions
Our analysis demonstrate that DEM based on Pléiades imagery are well-suited for the analysis of surface elevation changes on mountain glaciers, and especially offer the possibility for regional high spatio-temporal resolution assessments of remote areas for which no field data is available.As demonstrated, the relative accuracy is sufficient to investigate glacier surface elevation changes, even of magnitudes below 1 m, as experienced in the Khumbu Himal, and can therefore be applied over relatively short periods of time, such as those required for annual assessments.
However, some specific problems have to be mentioned and kept in mind while ordering, processing, and analyzing the data.Especially steep areas remain problematic: the noise levels are higher, and gaps occur more frequently, especially if the gradient of the terrain exceeds 40 • .For example, the DEM shows lots of gaps in the very steep south face of Nuptse or in the highest parts of Ngozumpa Glacier.The problem is compounded if steep areas are also affected by shadow, even though Pléiades offers a high image depth of 12 bit.This poses a particular problem for the study of geomorphological processes occurring on this steep terrain, such as the case of glacier accumulation via avalanches from steep headwalls.Fresh snow cover also poses a challenge for the DEM calculation, which is shown in the problems at the Ngozumpa glacier (see Section 2.2.1).
As no direct information on the mass can be gained from the remote sensing data, the conversion from volume to mass changes that is required for the calculation of geodetic glacier mass balances requires additional information.Especially for glaciers in areas where changing firn bodies can have a significant influence on the mass balance, as in the European Alps, or where the thickness and density of the firn bodies are unknown, as in the Khumbu Himal.Optical information can be used to determine the distribution of ice, snow, and less reliably, firn on the glacier surface during the time of acquisition, which can help with the determination of the spatial distribution of the increase or decrease of mass.However, ultimately assumptions, especially for the density of snow, are necessary to calculate glacier mass change from the volume change determined from the DEM differencing.Surface classifications of glaciers based on the high quality, multispectral Pleiades optical imagery can determine the debris cover or monitor snow and firn extent.This can be helpful to convert volume change to mass change, especially in areas where glaciers have changing firn bodies, and to extract surface objects, in order to map and monitor features such as proglacial lakes, ice faces, and avalanche cones that would help refine the geodetic mass balance calculations offered by the Pleiades workflow.

Figure 1 .
Figure 1.The study area in the Austrian Alps.One is the main part of Hochjochferner, which is investigated in this study.UTM 32 N.

Figure 1 .
Figure 1.The study area in the Austrian Alps.One is the main part of Hochjochferner, which is investigated in this study.UTM 32 N.
Remote Sens. 2017, 9, x FOR PEER REVIEW 4 of 26 sensed geodetic mass balances based on DEM differencing offers the best means of measuring current change in these glaciers.The role of the different processes contributing to the mass change of debris-covered Himalayan glaciers is still under investigation [30-32].

Figure 2 .
Figure 2. The study area in the Khumbu Himal; the red frame shows the acquisition area for the Pléiades data.The investigated glaciers are numbered, with 1 signifying the Ngozumpa glacier, 2 the Khumbu glacier, 3 the Nuptse glacier, 4 the Lhotse Nup glacier, and 5 the Lhotse glacier.Source: GeoEye via ESRI basemap.UTM 45 N.

Figure 2 .
Figure 2. The study area in the Khumbu Himal; the red frame shows the acquisition area for the Pléiades data.The investigated glaciers are numbered, with 1 signifying the Ngozumpa glacier, 2 the Khumbu glacier, 3 the Nuptse glacier, 4 the Lhotse Nup glacier, and 5 the Lhotse glacier.Source: GeoEye via ESRI basemap.UTM 45 N.

Figure 3 .
Figure 3. Overview of the autumn 2015 Pléaides images of the Khumbu Himal study area.UTM 45 N.

Figure 3 .
Figure 3. Overview of the autumn 2015 Pléaides images of the Khumbu Himal study area.UTM 45 N.

Figure 4 .
Figure 4. Overview of the spring 2016 Pléaides images of the Khumbu Himal study area.UTM 45 N.

Figure 5 .
Figure 5. Example of the three spring 2015 Pléaides image tiles of the Khumbu Himal study area.In the northern and eastern part, snow was still present during acquisition, which influenced the respective digital elevation model (DEM) processing.UTM 45 N.

Figure 4 .
Figure 4. Overview of the spring 2016 Pléaides images of the Khumbu Himal study area.UTM 45 N.

Figure 4 .
Figure 4. Overview of the spring 2016 Pléaides images of the Khumbu Himal study area.UTM 45 N.

Figure 5 .
Figure 5. Example of the three spring 2015 Pléaides image tiles of the Khumbu Himal study area.In the northern and eastern part, snow was still present during acquisition, which influenced the respective digital elevation model (DEM) processing.UTM 45 N.

Figure 5 .
Figure 5. Example of the three spring 2015 Pléaides image tiles of the Khumbu Himal study area.In the northern and eastern part, snow was still present during acquisition, which influenced the respective digital elevation model (DEM) processing.UTM 45 N.

Figure 6 .
Figure 6.Hillshade of the spring 2016 DEM showing part of the highest area of the Ngozumpa glacier.A problematic area due to snow cover is clearly visible in the area east of the ridge.

Figure 6 .
Figure 6.Hillshade of the spring 2016 DEM showing part of the highest area of the Ngozumpa glacier.A problematic area due to snow cover is clearly visible in the area east of the ridge.

Figure 7 .
Figure 7. Statistical function for 2014 Pléiades raw data and co-registered data of the Hochjochferner area.

Figure 7 .
Figure 7. Statistical function for 2014 Pléiades raw data and co-registered data of the Hochjochferner area.

Figure 8 .
Figure 8. Profile 1, located west of the tongue of Ngozumpa glacier.The profile shows the surface differences in off-glacier terrain between the spring 2015 and autumn 2015 acquisitions from south to north.

Figure 9 .
Figure 9. Profile 2, located east of the tongue of Khumbu glacier.The profile shows the surface differences in off-glacier terrain between the spring 2015 and autumn 2015 acquisitions from south to north.

Figure 8 .
Figure 8. Profile 1, located west of the tongue of Ngozumpa glacier.The profile shows the surface differences in off-glacier terrain between the spring 2015 and autumn 2015 acquisitions from south to north.

Figure 8 .
Figure 8. Profile 1, located west of the tongue of Ngozumpa glacier.The profile shows the surface differences in off-glacier terrain between the spring 2015 and autumn 2015 acquisitions from south to north.

Figure 9 .
Figure 9. Profile 2, located east of the tongue of Khumbu glacier.The profile shows the surface differences in off-glacier terrain between the spring 2015 and autumn 2015 acquisitions from south to north.

Figure 9 .
Figure 9. Profile 2, located east of the tongue of Khumbu glacier.The profile shows the surface differences in off-glacier terrain between the spring 2015 and autumn 2015 acquisitions from south to north.

Figure 10 .
Figure 10.Post co-registration surface difference and location of profile 1 and profile 2 in the area of overlapping scenes west of the tongue of Ngozumpa glacier and east of the tongue of Khumbu glacier.The green line shows the location of the profile, and the surface differences in areas outside the glaciers are clearly visible.The high differences in the areas of unfrozen lakes are clearly visible.

Figure 10 .
Figure 10.Post co-registration surface difference and location of profile 1 and profile 2 in the area of overlapping scenes west of the tongue of Ngozumpa glacier and east of the tongue of Khumbu glacier.The green line shows the location of the profile, and the surface differences in areas outside the glaciers are clearly visible.The high differences in the areas of unfrozen lakes are clearly visible.

Figure 12 .
Figure 12.One-year surface elevation change of the investigated part of Hochjochferner autumn 2014-autumn 2015.UTM 32 N.

Figure 12 .
Figure 12.One-year surface elevation change of the investigated part of Hochjochferner autumn 2014-autumn 2015.UTM 32 N.

Figure 13 .
Figure 13.One-year winter and summer surface elevation change of the Ngozumpa glacier.UTM 45 N.

Figure 13 .
Figure 13.One-year winter and summer surface elevation change of the Ngozumpa glacier.UTM 45 N.

Figure 17 .
Figure 17.Annual surface displacement on Ngozumpa glacier, normalized values based on the autumn 2015 and spring 2016 acquisitions.UTM 45 N.

26 Figure 18 .
Figure 18.Summer, winter, and annual surface elevation changes per 50-m elevation band of the five investigated glaciers in the Khumbu Himal.The highest parts of Ngozumpa glacier where data quality is problematic are marked in blue.

Figure 18 .
Figure 18.Summer, winter, and annual surface elevation changes per 50-m elevation band of the five investigated glaciers in the Khumbu Himal.The highest parts of Ngozumpa glacier where data quality is problematic are marked in blue.

Figure 19 .
Figure 19.One-year surface elevation change of the five investigated glaciers in the Khumbu Himal.These were calculated once from spring 2015 to spring 2016 and once combined with the surface elevation changes of the spring 2015-autumn 2015 and autumn 2015-spring 2016 periods.The highest parts of Ngozumpa glacier where data quality is problematic are marked in blue.

Figure 19 .
Figure 19.One-year surface elevation change of the five investigated glaciers in the Khumbu Himal.These were calculated once from spring 2015 to spring 2016 and once combined with the surface elevation changes of the spring 2015-autumn 2015 and autumn 2015-spring 2016 periods.The highest parts of Ngozumpa glacier where data quality is problematic are marked in blue.

Figure 20 .
Figure 20.Detailed view of the one-year surface elevation change at the highest parts of Lhotse glacier.Surface elevation change rates at the avalanche cones below the very steep north face of Nuptse are clearly distinguishable.UTM 45 N.

Figure 20 .
Figure 20.Detailed view of the one-year surface elevation change at the highest parts of Lhotse glacier.Surface elevation change rates at the avalanche cones below the very steep north face of Nuptse are clearly distinguishable.UTM 45

Table 2 .
Shifts of the slave DEM in relation to the master DEM calculated during the co-registration process.Multiple values for the Khumbu Himal study area are due to the three triplets needed to cover the whole study area.ALS: airborne laser scanning.

Table 3 .
Standard deviation of surface elevation difference over stable terrain.

Table 3 .
Standard deviation of surface elevation difference over stable terrain.

Table 3 .
Standard deviation of surface elevation difference over stable terrain.

Table 4 .
X and y offset of the ground control point (GCP) in relation to the 2016 DEM.

Table 5 .
Total volume change and average surface elevation change for the investigated glaciers and periods.

Table 5 .
Total volume change and average surface elevation change for the investigated glaciers and periods.