Reanalysis of multi-temporal aerial images of Storglaci¨aren, Sweden (1959–99) – Part 1: Determination of length, area, and volume changes

. Storglaci¨aren, located in the Kebnekaise massif in northern Sweden, has a long history of glaciological research. Early photo documentations date back to the late 19th century. Measurements of front position variations and distributed mass balance have been carried out since 1910 and 1945/46, respectively. In addition to these in-situ measurements, aerial photographs have been taken at decadal intervals since the beginning of the mass balance monitoring program and were used to produce topographic glacier maps. Inaccuracies in the maps were a challenge to early attempts to derive glacier volume changes and resulted in major differences when compared to the direct glaciological mass balances. In this study, we reanalyzed dia-positives of the original aerial photographs of 1959, -69, -80, -90 and -99 based on consistent photogrammetric processing. From the resulting digital elevation models and orthophotos, changes in length, area, and volume of Storglaci ¨ aren were computed between the survey years, including an assessment of related errors. Between 1959 and 1999, Storglaci¨aren lost an ice volume of 19 × 10 6 m 3 , which corresponds to a cumulative ice thickness loss of 5.69 m and a mean annual loss of 0.14 m. This ice loss resulted largely from a strong volume loss during the period 1959–80 and was partly compensated during the period 1980–99. As a consequence, the glacier shows a strong retreat in the 1960s, a slowing in the 1970s, and pseudo-stationary conditions in the 1980s and 1990s.


Introduction
Glacier volume change is of major concern when estimating the effects of climate change on factors such as sea level change and water resources in mountainous terrain (e.g., Lemke et al., 2007). With data from glacier change measurements spanning from century to decadal scales, techniques for measuring and processing data have evolved with developing technologies. Thus, it is important to carefully assess older data sets and homogenize these to remove any bias introduced by changing methods or processing. Homogenization of data series has become a standard for instrumental climate data series (e.g., Böhm et al., 2001) and is gaining importance as available glaciological data series increase in length (e.g., Thibert et al., 2008;Huss et al., 2009;Fischer, 2010).
Glaciological research at Storglaciären has a long tradition reaching back to the late 19th century, when early photo documentations were compiled in the Tarfala Valley . In 1910, the first detailed documentation with photographs and terrestrial photogrammetry was carried out (Holmlund, 1996) and continuous length change measurements were initiated. Systematic in-situ mass balance studies based on the direct glaciological method (cf. Østrem and Brugman, 1991) began in 1945/46 by the Swedish glaciologist Valter Schytt, when the Tarfala Research Station (Stockholm University) was constructed. Based on these continuous long-term observation series, a wealth of scientific studies has been produced, especially related to glacier mass balance measurements (e.g., Schytt, 1981;Holmlund and Jansson, 1999;Schneider and Jansson, 2004;334 T. Koblet et al.: Determination of length, area and volume changes of Storglaciären reconstructions (e.g., Holmlund, 1987;Raper et al., 1996;Linderholm et al., 2007), modelling (e.g., Hock, 1999;Albrecht et al., 2000;Hock and Holmgren, 2005), hydrology (e.g., Hock and Noetzli, 1997), and ice temperature (e.g., Pettersson et al., 2003Pettersson et al., , 2007. In addition to the field observations, aerial photographs have been taken in decadal intervals since the beginning of the mass balance monitoring program. The resulting data were used to produce topographic glacier maps for the years 1949, -59, -69, -80, and -90, as well as for early, qualitative volume change assessments as summarized in Holmlund (1996). A first quantitative computation of glacier volume changes was done by Holmlund (1987). The author superimposed the maps of 1949, -59, -69, and -80, calculated changes in thickness at contour lines intersections and contoured at 5 m intervals, finally obtaining volume changes by planimetering. A first comparison with the glaciological mass balances showed large differences that were assumed to come from errors in the datum level of the maps (Holmlund, 1987). Albrecht et al. (2000) repeated the calculation of volume changes, this time based on digitized contour maps of 1959, -69, -80, and -90, interpolated to a regular grid for comparison with the glaciological mass balances. Again, major differences were attributed to inaccuracies in the maps and an error in the datum level. Also, the mapped glacier outlines and related area changes vary between different authors (cf. Holmlund, 1987;Holmlund et al., 2005), probably due to different interpretations of the glacier margins, where determining snow and shadow conditions present challenges. To complicate matters further, each map has its own history of basic data and assumptions, methods, revisions, and related errors (Holmlund, 1996). In order to reduce the errors of the maps and the related digitization processes, it is necessary to go back to the original aerial images and process them with consistent standard techniques.
In this study, we reanalyze the original aerial photographs of 1959, -69, -80, -90, and -99, applying a consistent photogrammetric processing strategy for all survey years. The resulting digital elevation models (DEMs) and orthophotos are used to determine and analyse length, area, and volume changes of Storglaciären. The main objective of the present study is to find out whether a thorough reanalysis of existing aerial photographs allows for a more accurate quantification of glacier changes. Given the limited number of long volumetric mass balance series, the quality control of existing data is very important. Therefore, the results are discussed including a qualitative and quantitative assessment of related uncertainties, in view of results from this and earlier studies as well as resulting consequences for future research. In a second paper (Zemp et al., 2010) the photogrammetrically derived volume changes of the present paper are compared to the in-situ mass balance measurements.

Study site
Storglaciären (67 • 55 N, 18 • 35 E) is located in the Kebnekaise massif in northern Sweden. As described by Schytt (1959) and Østrem et al. (1973), Storglaciären is classified as a small valley glacier of about 3 km 2 (cf. Table 6) with an elevation range from 1130 to 1700 m a.s.l. and is characterized by a branched accumulation area (Fig. 1). Based on its thermal regime the glacier is described as a polythermal glacier with a perennially cold surface layer in the ablation area (Hooke et al., 1983;Holmlund and Eriksson, 1989;Pettersson et al., 2003). Based on radio-echo sounding maps of the bed topography and the 1959 surface map, the ice volume of Storglaciären is approximately 300×10 6 m 3 (Björnsson, 1981;Eriksson et al., 1993;Albrecht et al., 2000). Based on these data the average ice thickness amounts to 100 m, with maximum values of 250 m (Björnsson, 1981). The mean annual air temperature  at the Tarfala Research Station (1130 m a.s.l.) is −3.5 ± 0.9 • C; the average summer temperature  is 5.9 ± 1.2 • C with maximum temperatures up to 20-25 • C; the mean winter temperature  is −6.6 ± 1.1 • C with minimum temperatures of c. −25 • C (Grudd and Schneider, 1996, updated with unpublished data of Tarfala Research Station). The mean annual precipitation amounts to 1000 mm a −1 (Holmlund and Jansson, 2002; based on data from Tarfala Research Station since 1989) which must be considered as a minimum estimate with respect to Storglaciären.  -69, -80, -90, and -99 (Holmlund, 1996). Aside from the photographs of 1949, which were not available, the data from the years 1959-99 built the substantial database for the determination of length, area, and volume changes in this study. Since the original negatives of the aerial photographs were not available, copies of the positives were used for the photogrammetric analysis. Each survey consisted of three to four photographs with an overlap of 65% and therefore allowed for processing the images photogrammetrically. The mean scale of the aerial images is 1:30 000 (1:22 000 for the 1960 survey). An overview of the available data sets, including scales and survey dates, is given in Table 1.
For the photogrammetric DEM computation, a combination of analytical and digital photogrammetry was applied. The analytical work was performed on a WILD BC2000S plotter. The interior orientation was calculated from the fiducial marks given in the photographs and the camera calibration files, which were partly available at the Lantmäteriet. The 1999 survey was chosen as the master survey, since camera information and corresponding calibration files were identified clearly for 1999. This information was also available for 1990. Due to missing camera information on the earlier photographs, the corresponding parameters had to be assumed based on the available calibration files from Lantmäteriet.
The digital calculation of the DEMs was conducted within the photogrammetry software SocetSet (BAE Systems version 5.4.1). For this purpose the scanned photographs (resolution 12 µm) were used in combination with the interior orientation from the analytical step. To optimize the quality, an additional 53 tie points were identified for the relative orientation. Finally, the absolute orientation was performed by the implementation of 21 Ground Control Points (GCPs) from field surveys (Jansson and Pettersson, 1997, Tarfala Research Station, unpublished data). The pixel size of the resulting DEMs is 5 m. Based on these DEMs, orthophotos of each survey data were created within SocetSet.

Differential GPS survey (reference data)
In order to quantify the DEM accuracy by the comparison with independent reference data, 26 points were measured in the field with a TRIMBLE 4600LS differential global positioning system (dGPS). A first dGPS survey had been conducted already in the 1990s (Bomark and Lundberg, 1995;Jansson and Pettersson, 1997), but as only a few of the points are within reasonable distance of Storglaciären, additional Table 1. Orthophoto information. The survey dates are based on information from Lantmäteriet and are also labelled on the image margins of the original photographs for 1980/90/99. The mean scale of the images is calculated from the average flight altitude. The fraction of areas influenced by snow and shadow is determined within a bounding box around Storglaciären (cf. Fig. 2). The contrast at the glacier surface is qualitatively estimated from the orthophotos.  Bomark and Lundberg (1995) as well as Pettersson and  give comprehensive information on technical details and precision of the dGPS system used at the Tarfala Research Station. Fig. 2 shows the spatial distribution of the reference points. They are concentrated in the Tarfala Valley and on the mountain ridge on the southern side of Storglaciären since they had to be placed in non-moving terrain and ideally on bedrock.

Glacier mapping and calculation of changes in length, area and volume
Glacier mapping serves as an important basis for determining glacier geometry and its changes and is influenced by the available data as well as by the operator's knowledge and experience. In this study, mapping was conducted by stereoscopic interpretation of the orthophotos as well as by geomorphometric analysis of the DEMs (e.g., shaded reliefs, profiles). Visualizations and spatial calculations were performed using ESRI ArcGIS Desktop 9.3. The glacier outline was mapped and digitized for every survey date by the interpretation of the orthophotos and the DEMs. Based on experiences reported from previous studies (Albrecht et al., 2000;Holmlund et al., 2005) and discussions (Hock et al., 2008), it was assumed that the real changes of the areal extent of the accumulation area are smaller than the differences due to the interpreter. Hence, outlines of the glacier tongue were digitized for all the survey years, keeping the outline of the accumulation area consistent (based on the 1990 image). Length changes were calculated on the basis of the corresponding orthophoto (Fig. 3). Thereto, a band of stripes with 50 m distance was drawn parallel to the main flow direction of the glacier (i.e. W-E). Length change was  calculated as the average length from the intersection of the stripes with the glacier outlines (cf. Kappeler 2006). Based on the outlines of the different survey years, corresponding areas and area changes were calculated. The changes in ice thickness were quantified by subtracting the DEM at the earlier date from the DEM at the later date. Volume change was calculated by adding the thickness changes between two survey years within the outline at the time of the bigger extent and multiplied by the cell size (i.e. 25 m).

DEM and volume change accuracies
Probable sources of errors in the final DEMs are the quality of the raw data, technical information on the surveys (camera calibration files, dates), as well as ground control points. In addition, topography (steep slopes) and meteorological conditions (snow cover) influence the accuracy of a DEM. If a DEM represents the reality insufficiently, calculations based on DEMs are not reliable (Wood and Fisher, 1993;Temme et al., 2007). It is therefore essential that the accuracy, which is regarded as the difference between a recorded value and the true value, is specified (Jones, 1997). Many publications have examined the different methods to evaluate the accuracy of DEMs (e.g. Heuvelink, 1998;Kääb 2005;Fisher and Tate, 2006;Kraus, 2004;Kraus et al., 2006;Thibert et al., 2008). The range of methods varies from qualitative evaluation of shaded reliefs to complex error propagation mod- elling (e.g. Li et al., 2005). Thibert et al. (2008) present a sound and comprehensive error compilation approach for photogrammetric analysis of glacier volume changes. Such a detailed approach cannot be conducted within this study, as not all parameters are known for the early survey dates (due to missing camera calibration files) and because analogue and digital steps are combined in the photogrammetric processing.
In the following, we qualitatively evaluate the new DEMs and orthophotos and use three different approaches to quantitatively assess the uncertainty of the DEMs and derived elevation changes: (a) a comparison of the DEM values with an independent set of dGPS points in non-glacierized terrain and statistical significance testing of the glacier changes (signal) versus the elevation changes in non-glacierized terrain (noise), using (b) T-Test and (c) Monte-Carlo simulation approaches. The first two methods (a, b) are used to interactively estimate systematic and stochastic uncertainties of the glacier thickness and volume changes.

Qualitative evaluation of DEMs and orthophotos
The aerial photographs available for Storglaciären are well suited for the computation of DEMs and orthophotos, even though a lower accuracy must be expected in areas with steep terrain or low contrast (e.g., due to shadow, snow cover). A review of the different photographs (see Table 1) shows good contrast on the glacier surface in 1959, 1980, and 1999. An extract of the orthophoto of the year 1999 and the corresponding shaded relief are shown in Fig. 2. The hill-shading of the DEM shows a very plausible representation of the real terrain. Regional artefacts in the accumulation areas of Storglaciären and Isfallsglaciären correspond to areas with snow cover and a discoloration in the orthophoto. In 1959, the whole area is covered by a thin snow cover and a large part of the glacier is lying in shadow, which is again represented as artefacts in the hill-shading of the DEM. In 1969, some prominent artefacts show up in steep and shadowed terrain.
In 1990, most of the glacier area has a good contrast but a lot of snow is lying in concavities such as in front of the glacier tongue. This signal is reflected in local artefacts, where snow cover, shadow, and steep terrain occur in parallel. The best contrast is given in 1980, where most of the area is free of snow and only few shadows occur. However, analyzing the corresponding shaded relief, a striking rugged topography occurs in the area northwest of Storglaciären, indicating a general tilt in the DEM. In addition, a linear artefact is visible in the orographic right part of the accumulation area of the glacier, which seems to be a vertically small but systematic step in the DEM.

Comparison with dGPS reference points
Common methods to compare the elevation of a set of independent reference points with the elevation at the corresponding coordinates in a DEM are the calculation of (a) the mean value (MEAN), (b) the root mean square error (RMSE), (c) the standard deviation (STDV), and (d) the standard error (SE) (Kääb, 2005;Fisher and Tate, 2006). Obviously, the comparison of precise point data with pixels (in our case with side lengths of 5 m) imports certain inexactness in the analysis. All statistical values for each survey date are given in Table 2 The validation with dGPS data was used as an estimate of systematic and stochastic uncertainties for the glacier thickness change within an observation period (T0-T1) as follows: where n is the number of dGPS points, z the elevation, and SE the standard error of DEM and dGPS elevation differences. The resulting uncertainties are given in Table 3.

Elevation changes in glacierized vs. non-glacierized area
As a standard method, a Student's T-test was performed in order to test if the ice changes (signal) on the glacier are significantly different from the error (noise) in the non-glacierized area which is assumed to be stable between two surveys. Separate evaluations were conducted for the accumulation area, the ablation area, and the total area and give overall p-values below 0.001 so therefore can be considered to be significant. The only exception was a p-value of 0.36 for the total glacier area in the comparison of 1990-1999. This can be explained by the small thickness change that occurred in that period. From the analysis of elevation differences in non-glaciated terrain, the mean difference between two DEMs can be considered as the systematic uncertainty for the volume changes of the corresponding time period.
where n is the number of non-glacierized DEM grid cells. The related standard error provides a measure for the corresponding stochastic uncertainty: where STD is the standard deviation of the non-glacierized elevation differences of the two DEMs and n the number of DEM grid cells. Note that for the calculation of the latter (i.e., the standard error) the elevation differences of the grid cells in the sample have to be independent. Due to the photogrammetric auto-correlation this is not given when using all grid cells. We addressed this issue by reducing the sample of n grid cells by applying a nearest neighbor re-sampling in ESRI ArcGIS Desktop 9.3. Under the assumption that the auto-correlation of pixels with 100 m (or 20 pixels) distance is negligible, we used an n of 523. The resulting uncertainties are given in Table 3.

Slope dependency of errors
The accuracy in DEMs varies due to their dependence on good contrast. Possible problems during the photogrammetric process may occur in areas with steep slopes or strong relief (e.g. Gousie, 2005;Temme et al., 2007). A T-test as performed above simply compares glacierized and nonglacierized area and does not consider this slope dependency www.the-cryosphere.net/4/333/2010/ The Cryosphere, 4, 333-343, 2010 of the error. To account for this circumstance and to quantify the inaccuracies that accumulate in the results of numerical modelling, a Monte-Carlo Simulation was performed (Burrough and Mc-Donnell, 1998). It is assumed that the error in a DEM has a Gaussian (normal) probability distribution function with known mean µ and standard deviation σ . As the error is normally distributed µ is 0, σ varies through the DEM depending on the slope of the terrain and hence σ increases in steeper slopes. A linear relationship between the (absolute) differences in elevation of two DEMs in terrain where no changes are assumed (the glaciated areas are masked out) and the slope of the terrain can be calculated. Based on this relationship, σ is estimated performing a Maximum Likelihood Analysis. Summarized, the distribution of the error in a DEM is N 0,a · slope where a indicates the relationship between slope and error. Table 4 gives an overview of the used parameter for all time periods. The Monte-Carlo Simulation was then performed on the area of Storglaciären. With 10 000 iterations a random error was calculated for every pixel following the Gaussian probability distribution of the error (Eq. 5); following Burrough & McDonnell (1998) at least 100 iterations are required. In a next step the mean error was calculated for each iteration.Then the statistical comparison of the artificial data sets (noise) with the calculated thickness changes (signal) between two years showed the probability that the real changes on the glacier are equal from the random error in the DEMs. Table 4 specifies the coefficient between slope and error (a), as well as the probabilities for accumulation, ablation, and total glacier area of all time periods. The coefficient a shows that the dependence of the error on slope is higher in the periods including the DEM of 1980. The probability that the thickness change over the entire glacier equals the random error in the DEMs is 0.03 for the entire period of investigation, and can thus be rejected based on a significance level of 0.05. It is striking how the p-values for the periods 1969-1980 and 1980-1990 are considerably higher than in the remaining time steps, although the signal is relatively strong. Splitting the results of the Monte-Carlo Simulation into accumulation area and ablation area, the p-value is notably higher in the accumulation area. Thus the probability that thickness change is indistinguishable from random error is higher in the accumulation area. Otherwise the p-values in the ablation area are lower than calculated for the complete glacier area. Again it is clearly visible that the periods including the DEM of 1980 have higher p-values.

Accuracies used for change assessment
We used different statistical methods in order to assess the accuracy of the multi-temporal DEMs. In general it was shown that the signal-to-noise ratio (glacier changes vs. changes in non-glacierized terrain) is appropriate. Two different analyses were performed for the quantification of systematic and stochastic uncertainties: the comparison with independent dGPS measurements and the comparison of glacierized and non-glacierized regions in the data sets. The latter approach integrates all errors, such as general shifts, included in the DEMs, since no independent data set, such as a master DEM, is available. Therefore we decided to use the comparison with dGPS reference points, providing an independent validation, as appropriate measure to assess the quality of our data. Hence, the resulting values (Table 2) are implemented as integrative estimates in the comparison of volume changes with in-situ mass balance measurements (Zemp et al., 2010).

Length, area, and volume changes
Changes of Storglaciären in length, area, thickness and volume are given in Table 5. Variations at the glacier snout between 1959 and 1999 are illustrated in Fig. 3. During the first two decades the retreat of the glacier is clearly visible. After 1980 only small changes in glacier length occur and between 1990 and 1999 no length change could be measured. The retreat over the whole period of investigation amounts to about 180 m. The area changes between 1959 and 1999 amount to −0.11 km 2 or −3.3% (Table 5). Again most of the changes occur in the first two decades. After 1980, the area loss is only 0.01 km 2 whereas between 1959 and 1980 the loss is 0.10 km 2 . An exact determination of the position of the glacier front after 1980 is hampered by perennial snow patches (see Fig. 3).
During the first decade (1959-69) the thickness change is -4.57 m which results in a volume change of −15.22×10 6 m 3 . The thickness and volume loss continues in the following decade  although it is slightly smaller. The second half of the investigation period  is characterized by a general volume gain. Between 1980 and 1990 the increase of thickness is 1.51 m. During the last decade the increase of thickness amounts to 0.68 m. For the entire period results a volume loss of 18.95×10 6 m 3 . Hence, according to the available DEMs, the volume loss at Storglaciären between 1959 and 1999 originates from the beginning of the period of investigation  and particularly from the first decade.
The spatial distribution of thickness changes at Storglaciären is shown in Fig. 4, giving the pattern for the different time steps as well as for the entire period. The elevation changes are shown for an entire bounding box around the glacier in order to provide an idea of the signal to noise ratio. Note that the non-glacierized terrain is located in the upper right corner, as well as on the left and lower sides of the bounding box (Fig. 2). The intense volume loss between 1959 and 1969 is the result of a thickness reduction over the entire glacier area. The highest values of thickness loss are achieved at the snout of the glacier and at a small patch at the southern edge in the middle of the glacier. On large parts of the glacier tongue there is a moderate thickness loss. In the accumulation area the distribution of gain and loss looks rather random. During the following decade (1969-80) the dominant component is again thickness loss in the ablation area and again the highest values occur at the snout of the glacier and in the upper part of the ablation area whereas the accumulation zone experiences a thickness increase. In contrast to the first decade the accumulation area is dominated by thickness increase. A noticeable edge divides the ablation area from the accumulation area. The same edge appears again from 1980 to 1990. The overall increasing thickness between 1980 and 1990 results from the thickness gain in the upper part of the ablation area. During the last decade there is a trend to thickness gain on the entire glacier. Particularly in the lower part of the ablation area, the thickness increase is clearly visible. Comparing the single periods, it becomes apparent that the thickness increase that occurs in the accumulation zone between 1969 and 1980 seems to flow towards the snout of the glacier in the following two decades (Fig. 4).
The combined analysis of length and area changes, representing the changing shape of the entire glacier as a reaction to volumetric changes, provides an integrative impression of glacier reaction to a change in climatic forcing. The change from a distinct loss in thickness/volume of the first two decades (1959-80) to a gain in the third and fourth decades (1980-99) results in a delayed reaction of the glacier extent in the form of a slowing down of the retreat over the entire four decades. The regional thickness gain of the 1980s in the lower accumulation area is about to reach the glacier tongue in the fourth decade. This does not lead to a glacier re-advance because of continued negative mass balances after 1995, as apparent from in-situ measurements.

Discussion
During the photogrammetric processing of DEMs and orthophotos of the survey years, inaccuracies appeared due to missing raw data, as well as shadows and snow cover on the photographs. The 1980 DEM in particular shows a rugged topography and large elevation differences compared to DEMs of 1969 and1990 in the non-glacierized areas northwest of Storglaciären, even if the contrast in the original data is very good and the fractions of snow cover and shadow are lower than in the other images. A possible explanation is that the DEM is partly tilted, which is supported by the results from the comparison with the dGPS points and leads to the poor performance in the uncertainty assessment. The comparison with the independent dGPS points indicates that all the DEMs seem to be systematically too low by a few meters, which must be attributed to an elevation bias in the GCPs (used in the photogrammetric processing), as already noted by Bomark and Lundberg (1995). Apart from the DEM of 1980, the corresponding standard deviations and standard errors of below 1.8 m and 0.4 m, respectively, are quite reasonable for comparing elevations of point measurements with 5×5 m pixels. Taking this into account, and also considering that the original photographs as well as some of the basic photogrammetric parameters were not available, the quality of the DEMs is in general well suited for deriving glacier changes. The "lesson learnt" in this context is that there is a need for (a) a proper documentation of all available data sets and their acquisition (including survey dates, raw data and processing steps), (b) the maintenance of existing networks such as the ground control points around Storglaciären, and (c) the compilation of a new high resolution dataset that provides a reference DEM and allows for an improved quality assessment (Kääb, 2005).
The photogrammetrically derived glacier length changes are compared with the in-situ front variation measurements as available from the WGMS (2008 and earlier issues; Fig. 5). The in-situ observations of the glacier snout are based on tape measuring from four to five points across the snout to a single fixed point  and with the use of a theodolite with a distance meter from one single fixed point . The results from the in-situ observations and the photogrammetry agree well in the three periods between 1969 and 1999. In the first period , the cumulative annual in-situ measurements result in a retreat of 81 m, which is 24 m less than calculated from the orthophoto. The larger value of the remote sensing method comes from the measurement stripes on the northern half of the glacier tongue that intersect at a low angle with the glacier outlines. This opens the door to a larger discussion on the most appropriate method for determining heterogeneous front variations along complex glacier terminus geometries (cf. Kappeler, 2006), which cannot be followed here. As already mentioned earlier, the determination of the glacier area and corresponding changes depends to a certain degree on the interpretation and generalization of the interpreter. This is especially true when the relative changes are small, as in the case for Storglaciären. The comparison with glacier areas as mapped Table 6. Areal extent of Storglaciären as published in different studies.
Year Holmlund Holmlund this (1987) et al. (2005)  by Holmlund (1987) and as used for the calculation of specific mass balances by Holmlund et al. (2005) does confirm this (see Table 6). Whereas the absolute differences are relatively small, the outlines by Holmlund et al. (2005) show a stronger area loss between 1959 and 1969 followed by a growing glacier extent until 1990, which contradicts the present study. The differences come from the different interpretations of the glacier outlines in the accumulation area and -due to the small absolute values -only has a minor influence on the calculation of the mass balance (Zemp et al., 2010). However, with the continuation of glacier changes this issue must be addressed in future.
The computed volume changes are compared to the results from earlier studies as shown in Table 7. The only good agreements are in first decade  with the original result by Holmlund (1987) and in the second decade  with the corrected one by the same author. Note that the corrected volume changes by Holmlund (1987) are adjusted to better fit the glaciological mass balance measurements and, hence, a direct comparison might not be appropriate. The comparison with the volume changes 1959-69-80-90 as calculated by Albrecht et al. (2000) shows large differences: whereas they find only half the volume change from this study in the first decade, their results are twice and almost four times our values in the second and third decade, respectively. The different reference areas -as discussed abovecan only explain a small portion of these deviations. The major deviation is best explained by an error in the datum level (Morgner and Hock unpublished data) and the different methodologies applied. Hence, all studies based on data from Albrecht et al. (2000) need to be reanalysed and conclusion need to be revised. Holmlund (1987) and Albrecht et al. (2000) computed the volume changes from interpolating contour lines of the topographic glacier maps using both analogue and digital approaches. The use of these maps, instead of the original aerial photographs, is an intermediate production step that introduces additional errors. Additionally, all these maps were produced at different points in time by different operators using different methodologies (Holmlund, 1996). A sound quantification of these errors -e.g. by analyzing the changes in non-glacierized terrain -is not possible, as the corresponding data is not (digitally) available (Holmlund, 1996;Albrecht et al., 2000). The results of this study, based on dia-positives of the original aerial photographs using consistent photogrammetric processing by the same operator, allow for a sound uncertainty assessment. We therefore consider the volume changes presented here -including estimates for systematic and stochastic errors -as the most consistent dataset, which is now ready for comparison with the results from the direct glaciological mass balances (cf. Zemp et al., 2010).

Conclusions
The presented approach of reanalyzing dia-positives of original aerial photographs of 1959, -69, -80, -90, and -99 with standard photogrammetric techniques resulted in a complete and consistent dataset of DEMs and orthophotos of the Tarfala Valley. Based on this new dataset, changes of Storglaciären in length, area, and volume are computed for the time periods between these surveys. The glacier lost 15 and 11×10 6 m 3 from 1959-69 and 1969-80, respectively. In the following two decades  a partial regain of the lost ice volume of 5 and 2×10 6 m 3 was found. Over the entire period from 1959-99, Storglaciären lost an ice volume of 19×10 6 m 3 . Averaged over the glacier area, this corresponds to a total ice thickness loss of 5.7 m, or to a mean annual ice loss of 0.14 m. The glacier reacted to this volume change with a slowing of its retreat, finally to pseudostationary conditions in the last observation period .
The uncertainty assessment shows that elevations of all the DEMs are systematically too low by a few meters, but with standard errors of below one meter. Thereby, the DEM of 1980 performed worse than the other DEMs. Statistical comparisons of the glacier changes with the "noise" in nonglacierized terrain prove the general significance of the results, with ablation areas performing better than accumulation areas. Again, the results including the DEM of 1980 perform poorer than the others. From the first to the last observation period, the absolute signal decreases and challenges the basic dataset and methodology.
The resulting length changes fit well to the cumulative in-situ observations. Only in the first decade (1959-69) do larger differences occur, which are probably due to the different measurement approaches. The absolute changes in area are small, but the relative changes differ between this and earlier studies. This is attributed to different interpretations of the glacier margins, especially in regions with shadow and/or snow cover. The major differences from earlier studies, however, are in the resulting volume changes of the glacier. Although difficult to quantify, we see the major cause as the heterogeneous methodology of earlier studies, which derived volume changes indirectly from topographic glacier maps of various origins. The resulting volume changes from Table 7. Comparison of the volume change [10 6 m 3 ] at Storglaciären from the results of Holmlund (1987), Holmlund (1996) and Albrecht et al. (2000), and this study. The results by Holmlund (1996) and Albrecht et al. (2000) are probably based on the glacier areas described in Holmlund et al. (2005).
Time period Holmlund (1987) Holmlund (1987) Holmlund (1996) This study (original data) (adjusted data) and Albrecht (2000) 1949-59 +9. the present study -although not free of errors -are based on a consistent re-processing of the original material, come with a sound uncertainty assessment, and allow a comparison with the in-situ mass balance measurements (cf. Zemp et al., 2010). The resulting data promote the importance of aerial photographs in glacier research, especially for quantifying length, area, and volume changes as well as for crosschecking in-situ measurements. Beyond that, the available DEMs and orthophotos can be used for investigations on the surrounding glaciers and for geomorphological mapping purposes in the Tarfala Valley. Looking into the (scientific) future, we would like to stress the importance of such decadal flight campaigns. We highly recommend compiling a new, high-resolution and high-precision data set that serves as a master data set for a new accuracy assessment, including independent spatial data. In this context, the spatial distribution and visibility of the ground control points around Storglaciären have to be improved and the date of the survey should be as close as possible to the in-situ annual mass balance measurements. With such a data set, glacier changes since 1959 can be quantified with higher accuracy and reliability. In addition, more analyses using the multi-temporal DEMs in combination with reference data and glaciological field measurements should be performed at the test site Storglaciären. By the application of modern techniques such as airborne laser scanning, problems related to shadows and snow cover, influencing the final data product, can be reduced.