2.1. Study region

The study region covers the glaciers in the Alps that are located in Austria, France, Italy and Switzerland (Fig. 1). In all four countries a few small glaciers are excluded either due to adverse snow conditions or local clouds. These will be added to the inventory when more suitable scenes become available. According to the WGI, the Alps host about 5150 glaciers covering 2910 km² with sizes between 0.02 and 87 km² (Reference Haeberli, Boösch, Scherler, Østrem and WallénWGMS, 1989). Most glaciers are small (80% are <1 km²), while only a few (0.6%) are >10km² but these cover 22% of the total glacier area. The highest elevation is found at Mont Blanc (Monte Bianco; 4808m a.s.l.), and mean glacier elevation is ~2900∓300 m. With an assumed total glacier volume of about 100–130km³ in the Alps in the 1970s (Reference Haeberli and HoelzleHaeberli and Hoelzle, 1995; Reference Paul, Huggel and KääbPaul and others, 2004a), an average mean glacier thickness of about 35–45m can be derived. Of course, this average is strongly influenced by a few very large glaciers, and most glaciers (by number) have mean thickness values of 25 m or less (Reference Farinotti, Huss, Bauder and FunkFarinotti and others, 2009). Nearly all glaciers in the Alps are temperate, with the exception of glaciers in or originating at high altitudes (>3500 m) that can be partly cold (Reference Suter, Laternser, Haeberli, Frauenfelder and HoelzleSuter and others, 2001). Regarding precipitation conditions, most glaciers are situated in a transitional climatic regime, but many also experience more wet/maritime (at the northern rim) or dry/continental (in the interior) conditions. An overview of the locally highly variable precipitation conditions in the Alps is provided by Reference SchwarbSchwarb (2000).

Fig. 1. The European Alps study region as seen in natural colours with a mosaic of the Landsat 5 TM scenes used for the glacier mapping. Blue areas with black outlines indicate individual glaciers. Some regions with mostly small glaciers are not covered (e.g. Monte Visco, Zugspitze and Dolomites). The upper left inset shows the location of the study site in Europe (black square). The lower right image is a TM band 5, 4, 3 composite of the Silvretta Group (white square in main image) with clouds (white) partly hiding some glaciers (cyan). The yellow square denotes the location of Figure 2, the yellow circle the location of Figure 6a and the arrow the location of Figure 6b.

2.2. Datasets

In total we have analysed ten scenes from Landsat TM that were acquired between 30 July and 21 September 2003 (Table 1). All scenes were downloaded from the website http://glovis.usgs.gov where they were available at no cost after the opening of the US Geological Survey’s Landsat archive (USGS, 2008). We delayed the processing until all scenes were available with L1T geometric correction, i.e. orthorectified with a DEM. Scenes of the same path were mosaicked before the processing. In the Alps, the summer of 2003 witnessed exceptional cloud and snow conditions as two intense heatwaves occurred at the end of June and in the second week of August (e.g. Reference SchärSchär and others, 2004). Combined with the high temperatures in July, these caused all snow outside of glaciers and partly also on the glaciers to melt away, so we had near-perfect glacier mapping conditions in most regions during August. This is of particular value for product quality, as shadow regions are smaller than on the normally used late-September images due to the still high solar elevation in August. A snowfall event in early September caused poor snow conditions for the rest of the year, at least in most parts of the Alps. On some images a few glaciers were cloud-covered and we were fortunate to find suitable scenes (1 and 10 in Table 1) from neighbouring paths without clouds. However, for about ten glaciers in the Silvretta group (inset in Fig. 1) and Dauphiné we had to use additional scenes from 2009 to delineate their accumulation region.

Table 1. Overview of the ten used Landsat 5 TM scenes from 2003

Several smaller glaciers in various regions of Italy (e.g. Brenta group and Dolomites) are not yet included. To identify glacier entities and to make an overall calculation of area change we used the data from the WGI (Reference Haeberli, Boösch, Scherler, Østrem and WallénWGMS, 1989).

To derive topographic glacier inventory parameters and drainage divides, a consistent DEM for the entire Alps is needed (e.g. Reference PaulPaul and others, 2009b). Several empirical tests and direct comparisons of the SRTM DEM, the ASTER GDEM and the DEM25 from swisstopo for glaciers in the Swiss Alps (Reference Frey and PaulFrey and Paul, in press) resulted in the decision to use the SRTM DEM (the void-filled version) for this purpose. The main reasons are: its acquisition date (February 2000) was still quite close to the year of the satellite scenes and clearly defined; it has fewer artefacts over glaciers and more distinct ridges at correct elevations than the GDEM; and for the entire region the quality is more consistent. We downloaded the respective tiles from the Consultative Group for International Agriculture Research (CGIAR; http://srtm.csi.cgiar.org/) and re-projected the mosaic to Universal Transverse Mercator (UTM) zone 32 (World Geodetic System 1984 (WGS84) datum) with 60 m cell size (as a compromise between the 30m pixel size of the TM scenes and the 90m resolution of the SRTM DEM). At this resolution, the smoothing of the ridges is still acceptable and topographic inventory parameters remain nearly unchanged (Reference Frey and PaulFrey and Paul, in press).

A hillshade of this DEM is shown for a subregion of the study site in Figure 2a to visualize its level of detail.

Fig. 2. (a) Hillshade of the SRTM 3 DEM resampled to 60 m spatial resolution with glacier outlines in blue and drainage basins in black. (b) Colour-coded flow-direction grid with glacier outlines and drainage divides as in (a). Each colour represents one cardinal direction from north (pink), to east (blue), south (green) and west (orange).

The vector outlines for the glaciers of the Austrian Alps were recently made available on the University of Innsbruck website (http://imgi.uibk.ac.at/iceclim/glacierinventory). This dataset was manually digitized from high-resolution aerial photography mainly acquired between 1997 and 1999 (Reference Lambrecht and KuhnLambrecht and Kuhn, 2007). We downloaded this dataset and re-projected it to UTM zone 33 (with WGS84 datum) for overlay and comparison.

3.1. Glacier mapping

All glaciers are classified in the main processing stage with a thresholded ratio image (TM3/TM5) and an additional threshold in TM1 to improve the classification in shadow regions (e.g. Reference Paul and KääbPaul and Kääb, 2005). By dividing the high reflectance values of ice and snow in the visible band (TM3) by the very low reflectance values in the shortwave infrared (TM5), a strongly contrast-enhanced ratio image resulted that is used to separate snow and ice (i.e. glaciers) from other terrain with a threshold. Several studies proved that this simple method is also very robust (regarding the threshold selection) and accurate compared to other methods (e.g. PReference Paul, Huggel, Kääb, Kellenberger and Maischaul and others, 2003; Reference Andreassen, Paul, Kääb and HausbergAndreassen and others, 2008). Of course, glacier ice under optically thick debris cover is not mapped by this method, and misclassification of turbid lakes also occurs. Although semi-automated techniques to determine the debris-covered glacier parts have been developed (e.g. Reference Bishop, Bonk, Kamp and ShroderBishop and others, 2001; Reference Paul, Huggel and KääbPaul and others, 2004a), we have edited all outlines in the post-processing stage by manual digitization using contrast-enhanced false-colour composites of the respective satellite scenes in the background. In part, oblique perspective views in Google Earth^TM were used to facilitate the interpretation. Glacier outlines from all scenes were corrected in their original reference system (e.g. UTM33 for numbers 1 and 2 in Table 1) and only the final outlines were re-projected to UTM32 and mosaicked afterwards.

3.2. Calculation of glacier inventory data

As mentioned before, the DEM provides both drainage divides and topographic data for each glacier. For both purposes three additional raster datasets (grids) were derived from the DEM: flow-direction, slope and aspect grids. From the latter, sine and cosine grids were also calculated in order to determine a correct mean value for the circular variable aspect (e.g. Reference EvansEvans, 2006). The flow-direction grid was used to aid in digitizing a separate vector layer with the drainage divides of each glacier. In Figure 2b a colour-coded version of this flow-direction grid is shown for a subregion of the study site along with the glacier outlines and the drainage divides. For the Swiss Alps we used the basin vector layer that was already prepared for the Swiss glacier inventory 2000 (cf. Reference PaulPaul, 2007); for the glaciers in Austria, France and Italy we used the codes from the WGI (converted to a point vector file) to identify individual glaciers according to the previous inventory. This overlay revealed that several, mostly small, glaciers are not covered by our inventory, as they were either cloud-covered or have disappeared. Due to inconsistencies in the former coding of glacier entities (e.g. when different drainage divides were used), we have not tried to perform a glacier-by-glacier comparison with the WGI data.

The basin vector layer is then digitally intersected with the corrected outlines to obtain individual glacier polygons, and a country code is assigned to each glacier by digital intersection with the polygons for each country. These were manually adjusted in regions with glaciers to avoid spatial overlap. The topographic inventory data per glacier were calculated from zone statistics following the method described by Reference Paul, Kääb, Maisch, Kellenberger and HaeberliPaul and others (2002) and the guidelines by Reference PaulPaul and others (2009b) for each of the following parameters: minimum, maximum, mean and median elevation, mean slope and mean aspect. Glacier area and code is automatically provided by the GIS; length is not yet included but will be added later.

3.3. Accuracy assessment

We performed two types of accuracy assessment: (1) overlay with glacier outlines that were independently digitized for the new Austrian glacier inventory from aerial photography, and (2) by comparison with glacier extents as digitized from a high-resolution IKONOS image acquired 2 weeks after the satellite scene. For (1) we compared our outlines only qualitatively, as the snow conditions were too different for a quantitative assessment. However, we find this overlay very helpful as it provides insight into the quality of the automatically derived outlines and the interpretation of debris-covered glacier parts. For (2) we used screen shots from Google maps that were geocoded iteratively (assuming UTM projection, zone 32 and WGS84 datum) until a good match with the Landsat scene was achieved for the respective subset. The extent of three differently sized glaciers on the IKONOS image was digitized once and on the Landsat image five times (each by another analyst) to assess differences in interpretation. The mean value of the multiple digitizations and the area from the IKONOS scene is then compared to the automatically derived area.

4.1. Inventory parameters

The count and area per size class for all 3770 glaciers in the sample are depicted in Figure 3a in bar graph form using the relative values. Table 2 also lists the absolute numbers for all size classes. Though the size classes have somewhat arbitrary boundaries, they reflect well the general characteristic of most glacierized mountain ranges, with abundant small glaciers present in regions with a few larger ones. For this overall inventory we computed a total area of 2050 km², with 82% of all glaciers smaller than 0.5 km² (covering 15% of the total area) and 11% larger than 1 km² covering 75% of the area. The high number of glaciers (1030) in the 0.1–0.5km² class (27%) contribute a substantial amount (>11%) to the total area. Neglecting this fraction and the even smaller entities would thus result in an incomplete inventory. When analysing the size class distribution per country, the fractions are roughly similar, with the largest glaciers found in Switzerland and a comparably higher proportion of them in the Swiss sample. Compared to the area reported in the WGI (2909km²), the area loss is ~30%. As most of the loss took place after 1985 (Reference Paul, Kääb, Maisch, Kellenberger and HaeberliPaul and others, 2004b), a mean area loss of 2%a^–1 results over this time period.

Fig. 3. (a) Histogram of the count and area distribution per size class for the entire sample. (b) Same as (a) but for each aspect sector.

Table 2. Count and area covered per size class for the entire dataset

When analysing the number and area of glaciers per aspect sector (Fig. 3b), we found that 52.5% of all glaciers face to the northern quadrant (northwest, north, northeast) whereas only 25.8% face to the southern quadrant (southeast, south, southwest). The area covered by the glaciers facing in these directions is rather similar (57.2% and 25.7% respectively). A slightly higher proportion of glaciers face east than face west. The comparatively large area covered by glaciers with southeast aspect is due to the mean aspect of the largest glacier in the sample (Grosser Aletsch). The count of aspect sectors for individual DEM cells gives more similar amounts for all sectors. This is because most of the cells from several of the largest glaciers (e.g. Grosser Aletsch, Fiescher, Miage, Rhone, Pasterze) are directed south rather than north. At the DEM cell scale, the asymmetry of glacier aspects (Reference EvansEvans, 2006) is thus less pronounced.

The variability of mean elevation with glacier aspect is depicted for all individual glaciers in the scatter plot of Figure 4a along with mean values per aspect sector. Apart from the already mentioned higher number of glaciers facing north, a tendency for lower mean elevations of north-facing glaciers can be seen. In the mean over the eight aspect sectors, glaciers facing south or southwest have 200m higher mean elevations than north-facing glaciers. There is, however, also a correlation of mean elevation with precipitation amounts, and the regions along the northern Alpine rim receive more precipitation and have predominantly north-facing glaciers. This autocorrelation must be considered when interpreting the lower mean elevations in climatic terms.

Fig. 4. (a) Mean glacier elevation vs aspect for all individual glaciers and mean value per aspect sector (red dots). (b) Mean glacier slope vs glacier area for the entire sample.

A scatter plot showing mean slope (derived from the DEM cells) vs glacier size is shown in Figure 4b. The tendency of larger (valley) glaciers to have lower mean slopes and of smaller (mountain) glaciers to be steeper is well known. The more interesting point is the increasing scatter of values towards smaller glaciers. This implies that smaller glaciers in a mountain range such as the European Alps can have nearly any slope and can thus be situated in all kinds of terrain. The reaction of these small glaciers to climate change is thus highly variable (e.g. due to the increase of topographic shading in steep valley walls) and the smallest glaciers do not necessarily have to disappear once temperatures increase. The large scatter of mean slope values for glaciers of the same size (e.g. 5–30° for 10km² size) also hints at a high variability of mean thickness values and the importance of considering slope in the related calculations (e.g. Reference Haeberli and HoelzleHaeberli and Hoelzle, 1995).

When plotting minimum and maximum elevation vs glacier size (Fig. 5a), we found clearly different values only for glaciers larger than 1 km². Glaciers smaller than this but larger than 0.1 km² show some overlap and are restricted to the elevation range 2000–4000ma.s.l. Glaciers smaller than 0.1 km² occur in nearly all elevation bands, which confirms that their reaction to climate change is difficult to generalize. There is no dependence of mean elevation on glacier size. The hypsographic curves (area distribution with elevation) summarized for each of the four countries are illustrated in Figure 5b. Apart from the different magnitudes of glacierization in each country, all four have most of the ice in the 3000 m elevation band (~100 m) which is in good agreement with the mean overall elevation (Fig. 4a). The smallest elevation range is found for glaciers in Austria, where the highest mountain (Grossglockner) does not exceed 3800 m, whereas some glaciers in Italy and France originate at the highest mountain in the Alps, Mont Blanc with 4808m elevation. Comparing the approximately similar hypsographies for Austria and Italy, a generally 100–200m shift of the distribution to lower elevations is visible for glaciers in Austria.

Fig. 5. (a) Minimum (black) and maximum (red) elevation vs glacier size for the entire sample. (b) Area–elevation distribution summarized for each of the four countries (see legend).

4.2.

Accuracy

The overlay of glacier outlines from the recent Austrian inventory and 2003 (this study) reveals several interesting issues (Fig. 6). At first glance, there is an overall good agreement of the outlines (particularly in the unchanged accumulation region of the larger glaciers). However, on closer inspection differences appear:

Fig. 6. Examples from glacier outlines as mapped here for 2003 (yellow lines) and for the Austrian glacier inventory for the ~1997–99 period (black lines). The background shows false-colour composites with bands 5, 4, 3 (as red, green, blue) of the TM scene used (1 pixel = 30 m). (a) Olperer Ferner (O) in the Zillertal Alps; (b) Habachkees (H) in the Hohe Tauern range. The arrow points to a debris-covered glacier that was not properly corrected in the 2003 inventory.

1. 1. There is a slight shift of the 1998 outlines compared to the satellite scene. This indicates that the projection routine in the GIS software used has problems with the projection of the Austrian coordinate system

2. 2. Strong area losses can be seen in both images, particularly for smaller glaciers of which several have even disappeared. This confirms the previously found continuation of the strong glacier area loss in the Alps since 1998 (e.g. Reference Paul, Kääb and HaeberliPaul and others, 2007).

3. 3. The shape of the outlines from the Austrian inventory is often much more irregular and has locally extended parts. This indicates that perennial and perhaps also some seasonal snowfields were included in the glacier area. The obviously strong area loss seen at some glaciers might also be caused by disappearance of snow rather than glacier ice.

4. 4. For some debris-covered, dirty or shadowed glaciers, the 2003 outlines have underestimated the ‘real’ glacier extent. This indicates that a more careful interpretation during manual editing is required when working at the pixel level and that outlines from a former inventory are very helpful for glacier identification, even when snow conditions were less good and glacier extents have changed in the meantime.

Our major conclusion from this comparison is that it is not advisable to directly calculate glacier area changes from the two datasets as they were obviously derived differently. Apart from the required correction of the 2003 outlines in some places, the differences in interpretation of the glacier outlines dominate the visible change.

The outlines from the multiple digitizations for the quantitative accuracy assessment are depicted in Figure 7 with the IKONOS image in the background. The overlay reveals that the TM-derived outlines are generally within the variability of the manual digitizations, i.e. the automated classification with TM provides a suitable mean value. The good agreement is confirmed when comparing the TM-derived areas for the three glaciers to the mean value of the respective manual digitizations (Table 3): The differences are <0.01km² for all three glaciers. The relative differences are within ∓1% for three comparisons and between 4% and –11% for the other three. The –11% difference from the manual digitizations is likely related to difficulties in accurately digitizing this glacier manually, as the standard deviation of the area differences is 12%. Overall, the glacier is ~4% larger with TM than from the digitization on the high-resolution image.

Fig. 7. Overlay of the automatic classification with Landsat TM (black), a manual digitization using the high-resolution IKONOS scene (cyan) in the background (screenshot from Google maps), and the multiple manual digitizations of the outline based on the TM scene (all other colours). Both glaciers are located in South Tyrol, Italy, close to Hintereisferner: (a) Steinschlagferner, (b) Schwemser Ferner.

Table 3. Glacier areas as derived automatically from TM and by manual digitization on the TM and IKONOS image along with the relative differences from the size derived by TM. The standard deviation for the multiple digitizations is also given (Std dev.). The bottom row provides sums and differences for all three glaciers

5.1. Glacier mapping

Finding the correct threshold for glacier classification and removing gross errors (lakes, rivers) is comparably fast for an individual scene (<1 hour). However, it has to be mentioned that work on the 2003 inventory benefited greatly from ideal mapping conditions (there was generally no seasonal snow outside of glaciers) in the extremely hot and dry summer of 2003 which was devastating for most glaciers. For ~30 glaciers in several regions (e.g. Vanoise, Glarner Alps, Silvretta and Ankogel group) clouds hid part of the glacier perimeter and this was corrected by integrating additional Landsat scenes (mostly from September 2009). As these corrections were mainly required in the accumulation region, we consider the area changes due to the scene acquisition 6 years later to be small.

The new glacier inventory for the entire Alps presented here has the clear advantage of being compiled within a comparably short period of time (7 weeks) compared to the previous inventory available from the WGI (Reference Haeberli, Boösch, Scherler, Østrem and WallénWGMS, 1989). However, due to the reduced spatial resolution of the satellite data used (30 m) compared to the WGI that was derived from aerial photography, the satellite product might be considered less accurate, particularly for small and/or debris-covered glaciers. On the other hand, the false-colour composites with TM bands 5, 4 and 3 (as red, green, blue) often show the (clean) ice extent more clearly than panchromatic images and are much easier to interpret in this regard. New techniques like the determination of glacier extents from (repeat) airborne laser scanning (Reference Abermann, Fischer, Lambrecht and GeistAbermann and others, 2010) are promising, but still very expensive when applied over large regions.

Compared to the compilation of the previous alpine-wide glacier inventory from the mid-1970s, the creation of an inventory from automated classification of satellite data with manual editing and DEM fusion is certainly much faster. However, the workload to accurately digitize debris-covered glacier parts is high (up to two-thirds of the total time), and the work must be done repeatedly, i.e. for each inventory. The effort required to digitize the drainage divides to separate contiguous ice masses is also high (maybe one-third of the required time), but they need to be digitized only once. This vector layer can then also be used for the next inventory, given that the orthorectification of the satellite imagery is sufficiently accurate (e.g. Reference Svoboda and PaulSvoboda and Paul, 2009). For glaciers in Switzerland we have only reprojected the basins used for the 1998/99 inventory (Reference PaulPaul, 2007) to UTM zone 32 and reapplied them to the glacier map.

5.2. Accuracy considerations

The comparison with the last inventory for the Austrian Alps revealed general differences in the interpretation of glacier extents, making a quantitative calculation of area changes impossible. Though the consideration of perennial snowfields was required to be consistent with the interpretation of the 1969 inventory, the area changes derived by Reference Lambrecht and KuhnLambrecht and Kuhn (2007) are likely a lower bound as the included perennial snowfields do not change much over time (Reference DeBeer and SharpDeBeer and Sharp, 2009). If seasonal snow had also been included in some regions (as inferred from the shape of the outlines), the real area loss of the glaciers would be even higher.

In agreement with results from previous studies (e.g. Reference Paul and KääbPaul and Kääb, 2005; Reference Andreassen, Paul, Kääb and HausbergAndreassen and others, 2008), we found for the three glaciers analysed here in more detail a mean area difference of –1.5% compared to manually derived outlines on a high-resolution IKONOS image. For individual and, in particular, small or debris-covered glaciers the differences are likely higher, but in the mean of a larger sample they are likely not, as the deviations are normally distributed (i.e. omission and commisson errors average out as Fig. 7 shows). As the glaciers we analysed here are largely debris-free, we conservatively conclude that the accuracy of the automatically derived glacier area with TM is better than 3% in the mean for a larger sample of glaciers and up to 5–10% for individual glaciers. On the other hand, the standard deviations of the area differences from the multiple digitizations are also in the order of 5–10%, indicating that the automated mapping of clean ice is at least as precise as the manual digitization.

5.3. Glacier change

Although a direct comparison with the total area reported in the WGI has several shortcomings (e.g. different sample; reference date from the WGI spans a 30 year period), we consider the calculated area loss of 30% as a good approximation for the loss since the mid-1980s, as from the 1960s to 1980s most glaciers in the Alps showed little change in area or even advanced (Reference Paul, Kääb, Maisch, Kellenberger and HaeberliPaul and others, 2004b). It has also to be considered that we have not yet included several small glaciers in each country (e.g. around Monte Visco, Zugspitze and the Dolomites). We assume that their total area is in the order of 20–30km². Furthermore, we likely have failed to map the debris-covered parts of some smaller glaciers correctly (see Fig. 6b). We assume that this gives an additional 10–20km² of glacier area that is not included. Thus our total area is a lower bound, and a more correct value would be in the range 2080–2100km². The total area reduction since the mid-1970s is thus ~820km², or ~40km² (or 2%) a^–1 when 1984 is used as the starting year for the area loss. This is likely one of the highest mean annual area loss rates worldwide over this period.

5.4. Topographic parameters

The topographic parameters for each glacier were calculated from the SRTM DEM which was resampled to 60 m cell size during re-projection from geographic coordinates to UTM zone 32. Though this largely retains the original values in steep topography, a glacier with a size of 0.01 km² only covers about three DEM cells at this resolution. The derived elevation or slope/aspect values hence become somewhat arbitrary. Moreover, the quality of the SRTM DEM is variable, and the filled values that replace large data voids can be subject to considerable errors. We thus expect that for several glaciers the calculated topographic parameters are outside acceptable error bounds, particularly for glaciers smaller than 0.1 km². For glaciers larger than this, the derived values are likely sufficiently accurate, as a direct comparison of parameters as derived from two other DEMs has shown (Reference Frey and PaulFrey and Paul, in press). For parameters like minimum elevation, the close temporal coincidence of the DEM with the satellite image can actually be more important for its accuracy than the accuracy of the DEM.

The general east–west orientation of the Alps results in several north–south-oriented local valleys, with most glaciers (52%) having a mean orientation towards northern sectors. Thus, topography certainly plays a role or in some regions even constrains the dominant aspect of the glaciers. The 200m lower mean elevation for north-facing glaciers can also partly be explained by the higher precipitation values at the northern alpine rim and is thus not only a result of reduced global radiation and increased topographic shading. However, the correct consideration of these effects is difficult and they might be investigated in a future study.