The role of debris cover in the evolution of Zmuttgletscher , Switzerland , since the end of the Little Ice Age

Debris-covered glaciers often exhibit large, flat, slow-flowing tongues. Many of these glaciers show high thinning rates today despite thick debris cover. Due to lack of observations, most existing studies have neglected the dynamic interaction between debris cover and glacier evolution over longer time periods. The main aim of this study is to reveal this interaction by reconstructing changes of debris cover, glacier geometry, flow velocities, and surface features of Zmuttgletscher (Switzerland), based on historic maps, satellite images, aerial photographs, and field observations. We show that debris cover extent has 10 increased from ~13% to >32% of the total glacier surface since 1859 and that the debris is sufficiently thick to reduce ablation compared to bare ice over much of the ablation area. Despite the debris cover the volume loss of Zmuttgletscher is comparable to that of debris-free glaciers located in similar settings whereas changes in length and area have been small in comparison. Increased ice mass input in the 1970s and 1980s resulted in a temporary velocity increase, as well as a lowering of the upper margin of debris cover and exposed-ice area, and a reduction of ice cliffs. Since ~2001, the lowest ~1.5 km are stagnant despite 15 a slight increase in surface slope of the glacier tongue. We conclude that the debris cover governs the pattern of volume loss without changing its magnitude, which is due to the large ablation area and strong thinning in regions with thin debris further up-glacier and in the regions of meltwater channels and ice cliffs. At the same time rising temperatures lead to increasing debris cover and decreasing glacier dynamics, thereby slowing down length and area losses. 1 Motivation and Objectives 20 Debris-covered glaciers have been observed to show a delayed adjustment of their length to climatic changes (e.g. Ogilvie, 1904; Scherler et al., 2011). This behaviour can be explained by melt reduction due to insulation by the debris layer, which commonly increases in thickness towards the terminus (Östrem, 1959; Nakawo et al., 1986; Anderson and Anderson, 2018), and is expected to distinctly prolong the glacier’s response time (Jóhannesson et al., 1989). In some regions debris-covered glaciers exhibit similar rates of volume changes as debris-free glaciers (e.g. Kääb et al., 2012; Gardelle et al., 2012; Nuimura 25 et al., 2012). Proposed reasons for this behaviour are: the emergence of surface features with exposed ice (ice cliffs, water flow channels, and ponds), stronger thinning further up-glacier that compensates for the debris-induced melt reduction, reduced mass flux from the accumulation areas and decreasing emergence velocities at the tongue, and englacial ablation with subsequent collapsing of the glacier surface (Ragettli et al., 2016; Banerjee, 2017). Ice cliffs and ponds can enhance ablation in comparison to debris-covered surfaces and even debris-free ice and are common features on debris-covered glacier tongues 30 (Brun et al., 2016; Benn et al., 2012). Especially during periods of negative mass balance, the down-glacier increase of debriscover thickness reduces ablation through its insulating effect and can lead to a lower – and even reversed – mass balance gradient (Nakawo et al., 1999; Benn and Lehmkuhl, 2000; Benn et al., 2012; Ragettli et al., 2016; Rounce et al., 2018). Over time, the surface slope of the glacier tongue therefore often reduces, which can lead to a decrease in driving stress and thus ice flow velocity (Kääb, 2005; Bolch et al., 2008b; Quincey et al., 2009; Jouvet et al., 2011; Rowan et al., 2015). Furthermore, 35 with an increase in equilibrium line altitude (ELA) the englacial debris melts out faster, leading to an extended debris cover further up-glacier (Kirkbride and Deline, 2013; Carturan et al., 2013). As a result, strongly debris-covered glaciers often have long, flat, and low-lying tongues with low flow velocities or even stagnant parts at their tongues (Benn et al., 2012). Current The Cryosphere Discuss., https://doi.org/10.5194/tc-2018-292 Manuscript under review for journal The Cryosphere Discussion started: 18 January 2019 c © Author(s) 2019. CC BY 4.0 License.

They also came up with critical comments and major drawbacks, which mainly referred to the issue of losing the focus on the main message in the number of different analysis and the overall length of the manuscript. For the revision of the manuscript we thoroughly considered all the critiques and suggestions by the reviewers. Changes have been implemented in the revised manuscript, which has -as we think -improved substantially.
Following this overall comment, we first reply to the major concerns raised by the reviewers and present the changes that have been performed in the revised manuscript. They encompass a number of additional changes, i.e. text has been rephrased, added and removed (revised manuscript is ~2000 words shorter), restructuring of content and chapters, removal and re-organisation of figures. The changes performed are consistent with the ones we suggested in our first reply. Changes in reaction to specific comments are shown in the respective replies, as is our argumentation where we didn't follow the referee's advice or where we chose a different way of approaching the consideration. For further changes we want to refer the referees to the revised manuscript, which is provided as a track-changed word document and a pdf where all revisions are already included (see supplements zip file).

Reply to major concerns
R2: 'The title is a tad misleading'.
We agree with the concern of R2 that the title is not perfectly appropriate for the focus of the study. To better reflect the main content of the manuscript we chose to focus the title on the glacier evolution rather than on debris cover and adjusted the title to: 'Unravelling the evolution of Zmuttgletscher and its debris cover since the end of the Little Ice Age' R1 & R2: 'Text as well as figures are too long and partly repetitive. It is sometimes confusing to follow the story due to the large number of analysis and results'.
We acknowledge the fact that the manuscript was rather long in places and that the presentation of 17 figures significantly contributed to this length. We agree that the main messages of the manuscript got somewhat 'blurred' by the shear amount of data and analysis presented. We also agree that there would be enough information for two stories but the main point of the study is the longer-term glacier evolution and the interaction and feedback between the different variables. We therefore decided not to slice this paper as separating the results would not allow to bring these results together and investigate their inter-relationships (dynamic feedbacks).
Overall, we have made the text more concise, removed three figures completely and additionally three parts of figures, present some of the content in different format and even removed some smaller analysis which are not relevant for the main outcomes. Specific details of the performed changes are listed below: -Removal of the generation of a moraine-based DTM for the end of the LIA (part of chapter 3.5) as it was not used for further analysis. -Removal of the elevation change of different debris thickness classes (former chapters 3.6.2 and 4.4.3). The message of this analysis is the insulation effect of the debris cover, which is already contained in the analysis of ablation stakes. -Shift of Figure 3b (the table) to the supplements. It is not essential in the main text (but still useful to refer to in the discussion). -Removal of Figure 6. The debris surface is described in the text and the photos have been placed in the supplement. -Substantial shortening of the chapter on tongue-wide elevation changes (new chapter 4.6), which is now also easier to follow. -Modification of the concept and figure (Fig. 11) of the mass balance gradient (incorporated into chapter 4.6). As suggested by R2, we have done the analysis for elevation bins of 50 m instead of the glacier sections as in this way it is easier to understand the concept and to present the results. Due to this modification of the concept, do not anymore need the section map (former Fig. 10b) and the table (former Fig. 10c). We have reduced the number of periods from 10 to 7. Also the terminology has been changed from mass balance gradient to elevation change gradient, which is the correct term in this case. -The temperature time series (former Figure 14) has been integrated into Figure 3 (mass balance) as a small, additional subplot sharing the same x-axis. In this way the evolution of the geodetic mass balance can be directly compared to the temperature evolution and Figure  14 could be removed. -Former Figure 13 (flow velocities from feature tracking and radar interferometer) has been moved to the supplements. This figure displayed partly redundant information (reduction over time, stagnant frontal part; already contained in Fig. 12) and partly information that is not necessary for the main message. -Removal of Figure 14 (temperature time series, see comment above on integrating it in Fig.  3). -Strong reduction of the discussion on the comparison to other glaciers outside of the Alps.
Nevertheless, we have kept a basic comparison to give some useful context, e.g. very few studies have compared debris extent changes over time. -Complete re-organisation of the discussion. It is now structured in the main topics, as they also contain the main messages. The discussion has been shortened by approx. one fourth.
Overall, this resulted in the removal (or shift to supplement) of three figures (Fig.6, Fig.13, Fig.14) and three figure parts (3b, 10b, 10c) as well as the shortening of a significant amount of text (~2000 words).
R1 & R2: 'Text on analysis and results and even the discussion do not always focus on the main messages, making these messages a bit blurry.' First of all, we have made the main messages stick out clearer in the discussion and conclusion. To make it clear also here, these main messages are: -Climate is responsible for the glacier-wide mass balance -The observed debris increase over time is very strong and unprecedented -Debris cover provokes spatial and temporal change patterns (reduced length and area change) -Ice cliffs are persistent as features but do not compensate the debris insulation, thus the larger glacier area is responsible for the comparable mass balance -Ice dynamics directly influences ice cliff as well as debris cover evolution We have also amended the text in the introduction, results, and discussion to serve more arguments to arrive at these messages. To be clearer about the messages, we have introduced additional subchapters in the discussion that specifically contain the main messages. The additional sub-chapters also help to structure the discussion in a more logical way. Together with the language comments by the reviewers as well as an additional correction by a native-speaking glaciologist we are confident that the story is much clearer now, thus increasing the impact of the study.
R2: 'Methods are not always reproducible. ' We further clarified the methods and added additional details. For the technical aspects of the production and detailed analysis of the DTMs (specific comment) we refer to Mölg & Bolch (2017) as they described the methods in detail. Also note that this study focuses on the glacier evolution and the links between variables and is not designed as a technical assessment or development of methods.
Specific steps that have been undertaken regarding this concern are: -Slight changes in text by adding further details of explanation and amendments of the text for clearer language (e.g. for the geodetic mass balance, debris cover mapping, ice cliff influence) -Referring to further existing studies where references were not sufficient (e.g. in surface feature mapping) -Removal of analysis that are not necessary and thereby reducing the complexity of the overall study as well as specific text sections and analysis descriptions (e.g. the morainebased surface topography) -Change of analysis to a clearer concept (ice cliff influence) -Clarification of our terminology to better distinguish between exposed-ice areas and ice cliffs.
R2: 'Ice emergence and horizontal displacements should be considered when calculating the influence of ice cliffs on total elevation change.' The objective of this analysis is to understand the potential influence of ice cliffs on elevation change.
To do so, we compared area-averaged elevation change with or without ice-cliffs and their surrounding areas but could not see substantial differences. Importantly, we do not explicitly compare different time periods. Also, both horizontal displacements and cliff backwasting are inherently accounted for by applying the buffer. Thus, our main conclusions only refer to elevation change per date (not ablation) and are in our view robust and do not require emergence velocities and horizontal displacements. We have adjusted the description to make the objective of the analysis and the conclusion drawn from it clear.
R2: '…The elaborate comparisons of Zmutt with other glaciers could be removed/reduced. … Some of it would be suited for a review on the state of Swiss glaciers.' We do not quite agree with the referee. A major aspect of this paper is the difference/non-difference in evolution/dynamic behaviour between debris-covered and debris-free glaciers. Thus, we consider it crucial to set the evolution of Zmuttgletscher in context with other clean-ice and debris-covered glaciers in the region (on the basis of simple measures such as length, area, and elevation change patterns).
R2 raised a specific concern regarding the analysis represented by Figure 17. This analysis clearly yields that debris-free and debris-covered glaciers can be separated based on elevation change patterns, and that Zmuttgletscher can clearly be grouped to the latter. The y-axis of the figure shows only elevation bins and no absolute numbers because the absolute elevation was normalised for comparison. This is necessary to account for the different elevation ranges of the glaciers. To make this clear we have changed the axis caption from 'Elevation class' to 'Normalised elevation class'. R1: 'The influence of ice cliffs is shown in both a table and a figure.' We decided to keep Table 2 in the main text. It shows information that is complementary to the figure, such as ice-cliff area and area share per date. However, we tried to incorporate it more in the discussion and focus on the cliff area reduction as a consequence of increased dynamics. This is an important result to show as it is related to dynamic interactions (see main messages).

R1: 'I suggest placing [section 4.4.3] after you have discussed the changes in surface velocities.' …
We agree and have adjusted the structure of the manuscript: The surface flow velocities have been moved up and newly represent chapter 3.4 instead of 3.7 in the methods, and chapter 4.4 instead of 4.6 in the results. We follow this suggestion by R1 because the flow dynamics is somehow inherent in the surface features and the mass balance and even referred to in the text.
Unlike announced earlier, the chapter on the elevation change gradient has not been incorporated into the chapter of geodetic mass balances. Because we only considered surface elevation changes the formerly used term 'mass balance gradient' was not correct and we've change it to 'Elevation change gradient'. With this title and the content it is a logical part of the chapter on tongue-wide elevation changes.
The new structure of the manuscript is as follows:  5. Ice cliffs and related elevation change 4.6. Tongue-wide surface elevation change and patterns 5. Discussion and interpretation 5.1. Methods suitability and shortcomings 5.2. Chronological evolution and process interactions 5.3. Unprecedented debris cover increase 5.4. The influence of debris cover on glacier geometry 5.5. Climate-driven glacier mass balance and flow dynamics 5.6. The role of debris-related surface features 5.7. Dynamic interaction between flow velocities, debris cover, and ice cliffs 6. Conclusions R2: 'To have an accurate curve the melt rates have to be normalized by the clean ice melt rate at the same elevation. Does normalized in this figure mean that you had a stake at the same elevation on a clean part of ice, representatively close to the debris stake, that was used to normalize each point?' We have changed the depicted values from 'normalised melt' to 'melt per day', since the clean-ice stake used for normalisation was not placed in a comparable climate as the other stakes (which would not be possible). We present the results of the 7-weeks period to proof and roughly quantify the insulation effect of the debris and not to establish an individual 'Östrem-curve' for Zmuttgletscher, which we also do not need for our further analysis. We think that our conclusion of melt reduction can be drawn because (i) debris of the observed thickness does substantially reduce ablation, (ii) the ablation rather depends on debris thickness than on elevation (there is zero correlation between elevation and melt, but an R² of 0.95 when we apply a logarithmic function to the measurements above 0 cm thickness), (iii) the measurements for two stakes with 5 and 25 cm, respectively, yield almost the same melt rate even though the stakes are located in quite different elevations.

R1: 'Elevation change vs. debris thickness information is repetitive.'
We were not exactly sure what the referee referred to but we reduced the amount of information in the text as well as former Figure 6. Additionally, the analysis of elevation change on different debris thickness classes has been removed as the temporal coverage on debris thickness information is limited.
P4 L10 - Table 1 shows the topographic maps, satellite, and aerial images, so move " (Table 1)" to after the satellite images. You could consider breaking this into two separate blocks of data like "topographical maps, …, and satellite images (Table 1) in addition to various field observations and long-term temperature measurements." To make the distinction that Table 1 is related to these 3 products even more apparent. New P4 L6. Thank you, good suggestion. Implemented as suggested.
P5 L6 -First time DTM abbreviation is used it should be spelled out fully. Done.
P6 L15 -"This resulted in a glacier area information for each data" doesn't make sense. Please clarify. Also, is this different than the first line of this section stating that glacier area was measured since 1859? It appears to be repetitive. New P6 L1. This might be a misunderstanding: it said "…in a glacier area information for each date" (not datA). We've adjusted the sentence to "This time series of maps and orthophotos resulted in a glacier area value for each corresponding date since 1859" to make it clearer.
P6 L24 -until 1997, and an additional data point… New P6 L10. The clause has been removed, because the GLAMOS measurements for Zmuttgletscher are not anymore used due to their problems of data gaps and reference positions before 1946.
P6 L32-33 -A bit unclear: were the two historic maps with debris cover symbols manually digitized as well? If so consider, "Debris cover extent for the orthophotos and historic maps that contained a debris cover symbol were manually digitized (Figure 2a and b)." New P6 L19. Yes, all were georeferenced and digitised. Changed as follows: "Debris cover extent for the orthophotos and historic maps was manually digitized." The fact that they contained the debris cover symbol is already mentioned before and repetitive. P6 L33 -… Siegfriend map (1879) was verified using two photographs… New P6 L21. Done.
P6 L34 -This information was valuable …, which was the region of the strongest changes. Switching from past tense in first sentence to present tense in this last sentence. Suggest keeping the tenses the same for the paragraph to make it easier to read. New P6 L22. True. We've changed the second sentence to past tense as well.
P7 L5 -A bit confusing. First, you should reference Figure 1, since they show where these were provided. Second, what do you mean by "for setting the elevation change observations into context"? Do you mean you were doing this to compare to the elevation change estimates from aerial/satellite imagery? Or were you measuring the surface mass balance to be able to break the elevation change in the aerial/satellite imagery into the surface mass balance and the flux divergence? I assume this will become clearer in the results, but it should be clarified here as well. New P7 L1. It was done to have an idea of the melt, and to complement the elevation changes from DTMs with shorter-term field measurements.
We've changed the sentence to "Ablation was measured at seven points on the glacier tongue (Error! Reference source not found.) during summer 2017 to better understand the influence of debris on melt and to complement the longer-term information on elevation changes." P7 L12-14 -A bit unclear as well. Do you mean "Debris thickness was measured via manual excavation along several transects perpendicular to the glacier tongue."? This seems to be the case from the figure and results. New P7 L8. Exactly. We've changed the sentence and added additional information: "Debris thickness data were collected in the field by manual excavation along and in between three transects perpendicular to the glacier flow direction in September 2017 (…)." P7 L15-16 -"used to put observations into context" does not provide any information to the reader. What observations were you trying to put into context? What context? Additionally, stations had to be close to what? To each other? To the ablation stakes? Please clarify. New P7 L13. So many, so justified questions  We hope it is clearer now. We've also added additional information that had before only been presented in the discussion. New sentence: "The selected climate stations had to be as close to the study site as possible, lie at similar elevations, and cover a long period. Since no single station fulfilled all requirements, we used the stations in Sion (484 m, ~35 km distant) and Col du Grand St. Bernard (2472 m, ~40 km distant)." P8 L28 -which reaches up to … New P8 L13. Done.
P8 L40-41 -consider consolidating the two sentences "Independent of the method, 100 m elevations bins (starting from 2100 m) were used to get a representative value for elevation change calculations in order to reduce the susceptibility of elevation change calculations to outliers due to the incomplete coverage of some areas of the glacier." New P8 L20. Thank you for the suggestion! We've actually rearranged the whole paragraph since the explanation of the processing was maybe not quite clear: "In areas of data voids or artefacts in the DTM, especially in higher elevations, no surface elevation change values could be calculated. In order to reduce the sensitivity to data voids and outliers, 100 m elevation bins (starting from 2100 m) were used for the elevation change calculation. In elevation bins without DTM coverage, the average of three extrapolation methods was used as the final mass balance value. The following three methods were used to fill the data voids: (1) a linear relationship between elevation and elevation change (ref), which is based on the strong respective correlation up to R²=0.93; (2) using the same relationship but additionally setting the elevation change values to 0 when they become positive in the highest elevations; (3) using the mean elevation change value of the uppermost elevations that contain data and applying it to the glacierised area above." P9 L22 -Do you mean, "change in debris thickness may impact the changes in surface elevation"? As part of the consolidation process of the manuscript, this analysis was removed completely. The analysis was based on these measurements and the DTM difference maps of several periods. From the orthophotos and the elevation change maps it is pretty clear where thin and thick debris are located (very similar pattern today and in the past), as well as bare ice. But the uncertainty of the thickness classes for the earlier periods is large and therefore we decided to remove this analysis.
The message of this analysis is also contained in the ablation stake measurements, i.e. debris melt reduction.
P10 L18 -this study. Thank you for spotting. Done.
P10 L11 -The wording is a little awkward. Consider "From close to the end of the LIA (~1850s) to 2017, Zmuttgletscher has retreated by 1907 +/-12 m (12.1 +/-0.09 m/yr). The maximum rate of retreat peaked between 1961 and 1977 at 21.7 +/-0.04 m/yr." Figure 3 -is there a reason Figure 3b, which is a table, is listed as a figure? What is the cause of the difference between Zmuttgletscher (this study) and GLAMOS 2018? They appear to be substantially different, albeit showing similar trends. New P10 L2. Thank you, we changed the sentences as suggested. Former Figure 3 (new: Fig. 4) has been split up: the length change was slightly adapted (to distinguish the curves more easily) and kept in the main text as a separate figure. The two tables were moved to the supplement, where they still provide some useful information. The difference between the two ZG curves cannot be well explained (as stated in the discussion 5.1). We agree that it is confusing and not necessary, since we have the length change information from our own data. Therefore we decided to remove the GLAMOS length record for Zmuttgletscher, as stated in P6 L10. P10 L25 -"was with" does not make sense. "Until 2013" also sounds awkward. Consider something along the lines of "At the end of the LIA 2.8 +/-0.2km2 of Zmuttgletscher was debriscovered, which has increased to > 5 km2 in 2013. During this time, the total glacier area has …" New P10 L14. Thank you. We've changed the sentence accordingly: "At the end of the LIA 2.8±0.2 km² of Zmuttgletscher were debris-covered (~12.9% of the entire glacier area), which increased to 5.03±0.1 km² in 2013." P11 L4 -"all parts of the glacier" is fairly vague, perhaps in "all tributaries"? New P10 L19. We've changed the sentence to "Generally, the extent of debris has expanded up-glacier in TMG and SBG." P11 L5 -the extent of debris has expanded … Done. P11 L5 -you refer to the debris-covered extent getting closer to Dent d'Herens, but isn't SBG (which is included in since you state "both" at the beginning of the sentence) getting equally as close to the headwalls of Dent Blanche? New P10 L20. Well observed, that was missing. We've added it to the sentence: "In both areas the debris extent has expanded to above the ice fall into the former accumulation areas and starts now close to the contributing rock walls of Dent d'Hérens (TMG) and Dent Blanche (SBG)." P11 L4-10 -I think it's important to be specific when referring to properties of the debris cover. In this paragraph you are discussing the extent of debris cover, so this should be made clear. The "debris cover grew strongly" could refer to the thickness of debris cover or the extent of debris cover. Since both are discussed in this paper, it's important to always distinguish between the two. Yes, this is important to avoid confusion. We have revisited all respective references and corrected them accordingly. Thank you! P11 L8-10 -This sentence is very confusing. What does "after the exposure of the rock wall" mean? Please clarify. New P10 L23. We've changed the sentence and hope that it's clear now: "By 2010, the ice fall at STG had thinned sufficiently to disconnect and expose the rock wall beneath, from which a small rock fall detached between 2010 and 2013 (Fig 6). This debris mound covered an area of approx. 170x300 m and is now slowly transported down-glacier and spreading laterally" P11 L18-20: difficult to read. Consider "In such cases, the typical base layer consists of finegrained material (sand and even silt) lying directly on the ice, which is overlain by a few centimeters of pebbles. Above this typical base layer, there is no specific sorting." New P11 L5. Thank you for this suggestion. We have changed the sentences to: "Field investigations reveal a homogeneous debris cover in some regions with stone sizes in the top-most layer mainly between 10 and 30 cm in diameter, and a much more heterogeneous debris cover in other regions, with pebble-sized stones to boulders in the metre scaleError! Reference source not found.. The typical base layer of the debris consists of fine-grained material (sand and even silt), and is overlain by a few centimetres of pebbles." P11 L12 -You state these are over several periods, but the caption states these measurements are from 05/07 -22/08 2017, which is a single period in time. Please clarify. New P11 L16. True. We have corrected the text accordingly: "Ablation measurements from seven locations and over a seven weeks period in summer 2017…".
P12 L20 -what does "number" refer to? The number of ice cliffs? New P12 L20. We've changed the sentence to: "Even though the glacier has become more and more debriscovered, we did not find any clear trend in area and location of ice cliffs (Table 3)." P12 L21 -"lower parts of the glacier tongue until the terminus" is redundant. New P12 L21. We've changed the sentence to: "They are almost exclusively located in the lowest part of the glacier tongue (yellow area in Fig9)." P13 L2 -delete "only very" Done. P13 L10 -delete "only" Done. P13 L19-20 -A bit redundant saying in some periods, then stating the period, and then doing the same later. Additionally, I'm not sure "strong" is the best choice of words to describe elevation change. It'd be better to be specific: "most negative" or the "surface lowering was greatest", something along those lines. For example, "The elevation change was most negative near the terminus between 1946-1977. Since 2001 though, the elevation change over the tongue has become more heterogeneous." New P13 L13-25. We have rewritten the whole paragraph and included your suggestions on terminology. We hope it is now more concise, easier to read, and less repetitive. P13 L22 -"Between 2879 and 1946, the average surface elevation change was -0.63 m/yr and most pronounced at the terminus. Surface lowering increased to …". Note: This is a bit repetitive of L19-20. See comment above. Repetitive content was removed.
P13 L26 -the tongue's average surface elevation change was almost balanced at -0.11 m/yr. New P13 L17. Sentence changed to: "Between 1977 and 1983 the tongue's average surface elevation was almost stable at change rates of -0.11 ±0.29m yr-1." P14 L1 -use either "surface lowering" or "surface elevation change". Previously, you've used "surface elevation change", so be consistent. This is especially true because you're using positive and negative values for elevation change, so surface lowering with a positive value makes this difficult to read. See comment above. Paragraph rewritten and suggestions implemented.
P16 L4 -do you mean the ratio of thinning for thick debris compared to bare ice? As part of the consolidation process of the manuscript, this analysis was removed completely (see comment above).
P16 L14 -please clarify what "no clear hint of an increasing balancing effect" means? There is no clear trend in the mass balance gradient flattening towards the terminus? remove the statement? Because further down we say there is an effect of the debris cover… New P13 L26. We have re-written most of the subchapter, among others removed this statement, which had even been somewhat contradicted by the text below. Additionally, we agreed with R2 that when we want to assess the elevation change gradient it makes more sense to look at the changes of elevation-based classes than of spatial sections. Therefore, we have changed the analysis to look at elevation changes of 50 m classes. Additionally, we have reduced the number of periods from eleven to seven. Thus, the additional map as well as the table are not necessary anymore. The whole analysis should now be easier to understand and the figure easier to read. P16 L24 -again, this appears to be a ratio. Ratios are not very intuitive. For example, if both regions were positive, they could still have the same ratio. Why are you reporting and discussing the ratio and changes in ratios? Is it providing information about how "connected" the upper and lower portions of the glaciers are? Is it used to understanding the dynamical state that the glacier is in? Please add a sentence or two to clarify what the ratio is actually telling us about the glacier. See comment above.
Section 4.5 -Is this referring to the "glacier-wide" mass balance? If so change the heading. New P9 L5. It refers to glacier-wide mass balance indeed. We think this should be inherent in the word 'geodetic', as it needs to be done glacier-wide to be able to calculate the mass balance. We made it clear that it is done glacier-wide in the respective methods chapter.
P17 L2 -After 1988, the mass balance became more negative again, even while the lowest areas on the tongue had a stable or positive mass balance. New P9 L10. Changed.
P17 L4 -When did this negative trend intensify? Rearrange this sentence if it was in 2001. New P9 L11. It intensified after 1988. We've changed the sentence to "The negative trend continuously intensified and after 2001 …". P18 L4 -add the year that this was done again for reference, so that it is clear in the text. It is obvious this this is an independent measurement, so you can delete "as an independent measurement" New P12 L14. Done.
P19 L31-L36 -this reads as a methods section, i.e., the reason you chose the methods you did as opposed to a discussion section. Consider moving this to the methods, although I leave this choice up to you. The other ones "suitability" of methods seem to discuss other studies, shortcomings, etc., which is why I don't think the other parts of this subsection are out of place. New P15 L25. We have moved parts of this paragraph up to the methods section (P6 L10) and also shortened it, but kept the part about uncertainties in place. We hope that the division is clearer now. P19 L39 -"Nevertheless we think…" are you referring to the long-term trend or are you saying that despite potentially underestimating the area, you still think the trend is significant? New P15 L25. Yes exactly. We have rephrased it to: "For 1961, the debris cover extent might be slightly underestimated as the orthophoto does not fully cover the ablation area and a thin snow-cover was present in higher areas. This could partly explain the slight reduction in debris extent from 1946 (3.98±0.01 km²) to 1961 (3.79±0.01 km²)." P20 L13-5 -"similar behavior", perhaps similar trend? "pointing out" also sounds awkward here. Perhaps "Zmuttgletscher has shown a similar trend indicating that its glacier-wide mass balance has foremost been governed by climatic changes rather than a response to changes in debris cover." New P18 L27 We've changed the whole paragraph about the climatic evolution and incorporated this statement in a different way. The closest sentence is: "The evolution of glacier-wide mass balance of Zmuttgletscher since 1879 is on century and decadal time-scales comparable to trends on other debris-free glaciers in the Swiss Alps, which are closely related to variations in climate (Fig 3, refs)." P20 L31 -Be careful with "higher" elevation change rate as this could imply a positive elevation change, but here it seems to be referring to being more negative? Please clarify. We have changed the terminology in the whole paragraph. It should be clear now.
P21 L3-5 -A bit confused here as to how the debris-covered area was stable if the upper margin of the zone was moved down-glaciers? Are you suggesting that the glacier advanced and that's what kept the debris-covered area constant? Otherwise, if the upper margin, which I assume is referring to the debris/ice extent interface moved downglacier, how could the debris extent remain constant? New P16 L19 Yes, the way the sentence had been phrased was confusing. What happened was that only in the central part of the tongue the clean-ice areas were pushed down and widened. We have revised the sentence to make this clear: "At the same time, the upper boundary of debris cover in these areas migrated slightly down-glacier ( Fig. 12a), but because the debris cover continued to expand in other areas, the total debris-covered area was roughly stable during this period." P21 L8-10 -"again" typically does not follow the verb, e.g., became more negative again, or the debris-covered area strongly increased again. Thank you. We have adjusted the word order in several places throughout the document. P21 L23 -and also likely thicker debris New P19 L31 Done.
P25 L18-21 -In Section 5.2.1 you state that the magnitude of mass change is governed by changes in climate as opposed to changes in debris cover. However, here you are stating that the different magnitudes in mass balance between Zmuttgletscher compared to glaciers in the Himalaya is likely attributed to the differences in debris cover. While this could play a role in the differences between the two, how do you know that this difference is not caused by differences in climate forcing?
Correct, it could also be the case for the Himalayan glaciers that debris governs the pattern rather than the magnitude of mass changes. The comparison to glaciers in the Himalayas has been cut down considerably. Among others, this statement was removed.
P26 L6-14 -Can this be consolidated? While it is nice to place the changes in debris-covered area into a regional context, the point gets a bit lost in this long list of every study, region, and change.
We have now strongly cut down this section and at the same time included more information about the absolute debris cover from other studies. At the same time we think that the point has become clearer. Most references have been kept for interested readers.
P27 L7 -These observations prove a … New P21 L20 Changed to "These findings suggest a clear and direct …"

P1L1
The title is a tad misleading, maybe. Yes, debris cover is incorporated in the analysis and a major component, but the main message is the evolution of the glacier. The results even show that the overall mass evolution of the glacier is not strongly affected by debris? I'm thinking something in the line of "Post-LIA evolution of Zmuttgletscher and its debris cover". Fair remark. We changed the title to 'Unravelling the evolution of Zmuttgletscher and its debris cover since the end of the Little Ice Age' P1L6 Often -> generally Flat. Debris-covered glaciers are generally not flat. They are hummocky due to the spatially variable melt rates induced by the debris. I think you mean gently-sloped or of low gradient. New P1 L8 True. We've changed 'flat' to 'gently-sloping'.

P1L7
Today -> at present New P1 L8 Changed to: 'At present, many of these glaciers show high thinning rates despite thick debris cover.' P1L11 Increased from approx. 13% to more than 32%? Provide specific numbers. New P1 L13 We are of the opinion that the abstract should provide a rough picture of the results and at the same time be as concise as possible (also required by the journal). Therefore we kept the rough numbers as exact numbers would also need the uncertainty range, which makes the sentence longer and less easy to read.

P1L15
Maybe provide numbers for the area and cliff changes. ~2005; ~1.5 -> Again be specific with you numbers. New P1 L19. Thank you. We've removed the tilde as these numbers are as precise as they can be.

P1L20
Why not just call it introduction? New P1 L25 Isn't this just a question of personal preference? I just think 'motivation and objectives' is a bit more explicit. In any case it takes up one line.

P1L25
Similar rates of thinning at the same elevation bands. Volume changes would imply mass balance but these are largely unknown from these studies. Also no reference to the (more recent) elevation difference paper here [Brun et al., 2017]? New P1 L30 Changed to "… exhibit similar thinning rates as debris-free glaciers…".
As far as we can see, Fanny Brun et al. 2017 did not distinguish between debris-free and debris-covered glaciers and didn't make any respective statements. Therefore it would not be correct to cite this paper in this context.

P1L26-29
There is some debate the last years about the importance of supraglacial ponds, ice cliffs and glacier dynamics/emergence. I think there are a number of recent papers that could be added here that touch upon this topic, e.g. [Pellicciotti et al., 2015;Vincent et al., 2016;Brun et al., 2018;Miles et al., 2018]. New P1 L35 Correct, we've added a few of these. There would be even more interesting studies, so we've also added an 'e.g.' before the references.

P1L37
Wouldn't call this 'as a result'. These glaciers also often have a long flat (low-gradient!) tongue because the insulating surface debris layer just allows them to extend into the lower, warmer valleys. New P1 L37 We have re-written the sentences so that it becomes clear that all these things are linked. The fact about the lower elevation of the glacier tongue is also linked to the reduced ablation, which was mentioned two sentences before we tried to summarise the above in the sentence that starts with 'as a result'. We hope this is clear now. We have changed it to: "Especially during periods of negative mass balance, the down-glacier increase in debris-cover thickness reduces ablation through its insulating effect and can lead to a lower -and even reversed -mass balance gradient (refs). Over time, this reduction in ablation can lead to a decrease in surface slope and, consequently, driving stress and ice flow velocity (refs). Furthermore, with an increase in equilibrium line altitude (ELA) the englacial debris melts out earlier, leading to an extended debris cover further up-glacier (refs). As a consequence of reduced ablation and driving stress, heavily debris-covered glaciers often have long, gently-sloped, and low-lying tongues with low flow velocities or even stagnant parts." P2L1-6 Long. Rephrase. New P2 L6-11 We have rephrased the sentence to: "Current research at many debris-covered glaciers and specifically in the Himalayas is mostly focussed on processes, such as ablation beneath the debris cover and in areas of ice cliffs and ponds, thinning of glacier tongues, surface flow velocities, or changes in debris cover over time (refs). A few studies have started to integrate the existing understanding and interactions of some of these processes into numerical models for the evolution of such glaciers (refs)." P2L7-12 Same here. Basically two paragraphs of one sentence each. Also these long itemizations could use some inline numbering. New P2 L12 Very good idea, we have added a numbering of the different shortcomings to distinguish between the main items of this paragraph: "However, most studies … of debris-covered glaciers: (i) time series are often short …; (ii) investigations are often local …; (iii) by considering only one or few variables … ; (iv) studies at glacier scale … ." However, grammatically, the paragraph consists of five 'sentences', not one. A semicolon can be used to separate independent clauses, the content of which is closely related. We consider this an element of style and less important for the understanding of the paragraph and therefore prefer to keep it.

P2L22
I think the importance part can be skipped, not quite relevant here and statement has an endless list of footnotes. Yes, true. We've removed it (" -which are in the center of attention due to their importance -").

P2L36
Are the latlon for the peaks relevant? The only important one, for the glacier, is missing. New P3 L2-5 The one for the glacier is in Figure 1a, but we've added it to the text. It's also true that the important information is rather the elevation of the peaks than their coordinates. We've removed them, making the text also a little bit nicer to read.

P2L40
Only originates from the rock walls? What about the other possible sources, see e.g. [Evatt et al., 2015]? P3 L5 You're absolutely right. In the discussion we also state that there are potential other sources (like moraines). We've added 'mainly': "The debris mainly originates from the surrounding rock walls." However, we've not added Evatt et al. (2015) because it's rather a theoretical consideration about the energy balance within the debris layer than actual analysis of debris sources. Therefore, in the discussion we stick to van Woerkom et al. (2018). P3L9-10 "at near-by almost"???
Can similar values really be assumed. We often considerably different behaviour of glaciers in a valley due to differences in microtopography etc. New P3 L8 We have changed the sentence to: "There are no direct measurements at higher elevations, but model estimates suggest values between 800-1500 mm (ref). Glaciological mass balance measurements at the nearby debris-free Findelgletscher (15 km distant) suggest end-of-winter accumulation around 0.8-1.5 m water equivalent (w.e.) (ref). However, at Zmuttgletscher, avalanching additionally contributes to accumulation." The measurements from Findelgletscher are the best guess we have at our disposition and we think it is valid to provide this information. However, we have rephrased the sentences to make clear that this is just a best guess and not the exact numbers we assume for Zmuttgletscher.

P3
I think that there is some irrelevant information provided in the study area section. Keep it clean and simple. We have taken out information about climate, geology, and other literature. Thus, we could shorten the chapter from 662 words to 452 words.
P5 section 3.1 Because the results depend strongly on the DTMs I think this section is a bit short. There are all kinds of caveats and things that can go wrong with DTM generation (especially working without fiducials and markings in Agisoft). A bit more detail on the exact procedures followed would be welcome. New P5, section 3.1 We have provided some further details. However, the generation of the DTMs is described in detail in a paper by Mölg and Bolch (2017), including a quality comparison of the output of two SfM software packages and a photogrammetry software (ERDAS Imagine). Among others, the idea of that paper was to have a base to refer to in order to save spae in the methods description in the current study.

P5L25
The DTM was produced FROM the tri-stereo image using photogrammetry, I'd guess.
New P5 L21 True, here the egg laid a hen. We've changed the sentence to: "In addition, we generated a glacier-wide DTM for 2017 (spatial resolution 1 m) using Pléiades tri-stereo, high-resolution satellite imagery based on photogrammetric principles using PCI Geomatica (ref)." P6L8 manually digitizing / manual digitization New P5 L29 Done.

P6L15
I think calling this ice fall is quite confusing terminology. Maybe ice deposits? New P6 L5 We've changed the sentence to: "…since they contribute mass to the main glacier through frequent ice avalanches." P6L17 …allowed correct interpretation of… New P6 L3 Done.

P6L29-30
Isn't this just normal error propagation? If it is not, why did you choose to New P6 L16. Yes, indeed. Isn't this what we say? The sentence is: "For the calculation of the changes, the uncertainties of the two respective dates are combined by error propagation, analogue to the mass balance data (ref).

P6L32
which -> that New P6 L19 We've changed the sentence to: "Debris cover extent for the orthophotos and historic maps was manually digitized (Fig 2)." We've removed the 'debris-cover symbol' since it is already stated before.

P7L5
Undertook -> performed New P7 L1 We've changed the sentence to: "Ablation was measured at seven points on the glacier tongue during summer 2017 to better understand the influence of debris on melt and to complement the longer-term information on elevation changes." for setting -> to put see above P7L6 …two metre long PVF stakes… New P7 L3-4 We've changed the word order. Why 'PVF'? I think the abbreviation 'PVC' comes from 'polyvinyl chloride'. But maybe this was just a typing mistake.

P7L12-14
When and how often were these measurements performed. Need more details. New P7 L8-11 The data was taken during a campaign in September 2017. We have changed the paragraph as follows and hope it is clear now: "Debris thickness data were collected in the field by manual excavation along and in between three transects perpendicular to the glacier flow direction in September 2017 (for transect locations see Error! Reference source not found.; results of upper transect see Error! Reference source not found.a). In general, each data point represents the average of three measurements ~1 m apart, and the standard deviation is used as uncertainty measure. For debris thicknesses above 20 cm only one measurement was taken."

P8L1-6
Have you not considered the object-based ice cliff mapping I've done before ? New P7 L39 Yes, we have. We have now included a reference to your paper in the methods section. The problem of a more sophisticated approach, as you described it, is that we had very variable input data, mostly black & white, with different texture, quality, and resolution. Also, the link to topography was very weak, i.e. there are ice cliffs with a relatively small slope (<25°) and steep areas (>40°) that are still debris-covered. The most important issues are discussed in the methods discussion chapter.
P8L17-18 I do not understand what you mean here. How can one assume a plane due to curvature? Or is the plane curved to mimic the laterally slightly convex glacier surface? But then it is not a plane, right? During the consolidation process of the manuscript this section has been removed because it had not been used for further analysis.

P8L32
I'm not sure if "stand" is the right word here Included in the removed sections, see comment above.

P8L36-40
It is completely unclear to me when each specific method was used and for what reasons. New P8 L19-29 It is simpler than it sounded: all methods were used for each mass balance period. For elevation bins with no coverage, the three extrapolation methods were applied and the average of the three methods was taken as final value, while the standard deviation went into the uncertainty calculation. We've changed the structure of the whole paragraph and hope that it's clear now: "Glacier-wide geodetic mass balance was calculated between dates for which DTMs covered large parts of the glacier surface (1879,1946,1977,1988,2001,2010,2017). In areas of data voids or artefacts in the DTM, especially in higher elevations, no surface elevation change values could be calculated. In order to reduce the sensitivity to data voids and outliers, 100 m elevation bins (starting from 2100 m) were used for the elevation change calculation. In elevation bins without DTM coverage, the average of three extrapolation methods was used as the final mass balance value. The following three methods were used to fill the data voids: (1) a linear relationship between elevation and elevation change (refs), which is based on the strong respective correlation up to R²=0.93; (2) using the same relationship but additionally setting the elevation change values to 0 when they become positive in the highest elevations; (3) using the mean elevation change value of the uppermost elevations that contain data and applying it to the glacierised area above. The area-weighted average elevation change of all bins was summed up to reach the average elevation change of the total glacier surface." P9L2-3 Just a complicated way to say area-weighted? New P8 L27 Correct. We've changed it to: "The area-weighted average elevation change of all bins was summed up to reach the average elevation change of the total glacier surface." P9L10-11 I do not understand why this was done. Now you're assuming a single elevation for the entire period for specific parts of the glacier? Wouldn't it be better to leave it as no data and perform weighted statistics appropriately?
We have removed this section since it adds complication and confuses and is irrelevant for the total uncertainty. For the record and in case you are interested: We did not use the elevation of 2010 for the entire period, which would indeed introduce a relatively big error. Instead, we only used the AREA of the elevation bins from 2010 and applied them to those dates, where the DTMs did not cover all elevation bins.

P9L22
impact on surface -> impact surface This analysis has been removed in the revised manuscript.

P9L23
So there is a class for 5 cm debris and for debris thicker than 15 cm. What about the rest? That is, between 5 and 15? Or do you mean just two classes, thin and thick, <15 and >15? As part of the consolidation process of the manuscript, this analysis was removed completely (see comment below).

P24-P26
So now there is suddenly a lot of debris thickness data. I don't understand this? How were the maps produced? Were there that many pits dug for all these time steps? If so, that's quite impressive. The only debris thickness measurements are from 2017. The analysis was based on these measurements and the DTM difference maps of several periods. From the orthophotos and the elevation change maps it is pretty clear where thin and thick debris are located (very similar pattern today and in the past), as well as bare ice. But the uncertainty of the thickness classes for the earlier periods is relatively large and therefore we decided to remove this analysis. The message of this analysis is also contained in the ablation stake measurements, i.e. debris melt reduction.

P9L35
Correlation quality ('strength')? I'd just use 'correlation'. Correlation itself implies a quality/strength of fit. Yes, absolutely. The term 'Strength' referred to the output parameter of this name by the software. It is confusing, thus we have removed it.

P10L8
Why smaller or equal than 0.03, and not just a uncertainty value? New P7 L36 Changed to ±0.03 m/yr.

P10L8
There should be no space between the first number and the plusminus symbol. You use m/yr semi-consistently throughout the manuscript. Preferebly use scientific notation and in glaciology it is somehow common to use pro annum instead of per year: 12.1 m/yr -> 12.1 m a-1. Change this throughout. We've removed the space throughout the manuscript. We've changed m/yr to m yr-1 (as it is used by many others, e.g. Benn et al. (2012), Fischer et al. (2015), Brun et al. (2017) etc.). This way the abbreviated words are consistently in English (which I find somehow more intuitive to read). P10L14 slowed down -> decelerated New P10 L5 Done. Figure 3 is a bit confusing with all the length changes and additional table. What is the point here. Maybe combine figure 3 and 4, and skip the display of the other glaciers? Just mentioning in text that they are quite similar around Switzerland would suffice, I think. New P10 L8 It is true that the combination with the tables was a bit confusing. We have removed the tables and shifted them to the supplements, where they still provide useful information (esp. on area changes). Nevertheless, we have kept the length change figure. It was slightly adapted to be easier to distinguish the curves. This length change curves of the different glaciers -which are actually quite different -are imortant for the interpretation of a debris-cover effect in the length change. Among others, this figure is used in the discussion to arrive at this conclusion. P10L27 …resulting in 33% debris cover at present. New P10 L15 Sentence changed to "During this time the total glacier area has decreased from 21.24±0.4 km² to 15.82±0.3 km², resulting in 31.8±2% of the glacier area being debris-covered in 2013 (fig)."

P11 Fig5
The different periods are a bit difficult to read with the current colour scheme. This could be improved by mapping it over a wider range of luminosities, i.e. from lighter yellows to darker reds. New P11 L1 We agree that the colour palette was not ideal. We've changed it to run from dark brown to orange to yellow.
P11L16-17 I'm not a geomorphologist, but is superficial the right terminology? "in the metre scale" sounds strange New P11 L5-7 Yes. We've changed the sentence to: 'Field investigations reveal a homogeneous debris cover in some regions with stone sizes in the top-most layer mainly between 10 and 30 cm in diameter, and a much more heterogeneous debris cover in other regions, with pebble-sized stones to boulders in the metre scale …'

P12 Fig7 panel b
No error bars on these points? The melt rate below the debris surface depends both on elevation (i.e. temperature + radiation) and on debris thickness. To have an accurate curve the melt rates have to be normalized by the clean ice melt rate at the same elevation. Does normalized in this figure mean that you had a stake at the same elevation on a clean part of ice, representatively close to the debris stake, that was used to normalize each point? If so, did you clear the glacier of debris or did you find a naturally clear spot? Could there then have been errors in the clean ice ablation measurements by radiation emitted or reflected by nearby debris; an energy balance component that would not exist on a completely debris-free glacier? EDIT: I see now in the text there was one reference stake at 2600. I'm then not quite sure whether the conclusion made are sound. New P11 L11 Correct, the reference stake used for the normalisation is the debris-free stake on 2600m. This is the lowest elevation with clean ice. Our conclusions are that (i) debris of the observed thickness does substantially reduce ablation, and that (ii) the ablation rather depends on debris thickness than on elevation. Looking at the individual stakes we see zero correlation between elevation and melt, but an R2 of 0.95 when we apply a logarithmic function to the measurements above 0 cm thickness. Also, the measurements for stakes with 5 and 25 cm yield almost the same melt rate even though the stakes are in quite different elevations. All measurements were taken on roughly flat surfaces. Since we do not use the derived Östrem-curve for further analysis, the exact numbers do not have any impact. But we think there is sufficient evidence to arrive at the conclusions mentioned above and in the text. We have changed to figure to show melt rates over the 7-weeks period. Thus, we arrive at the same conclusion without introducing a small additional error by normalising with the clean-ice stake from a different elevation. We also added a sentence about the uncertainties. They are too small to display as error bars in the figure.
P12L16-17 Could perturbation of the debris layer during drilling of the stakes have caused a difference in the debris matrix that could have affected the melt rates? We gave our best to replace the debris the same way it was before, but some deviation from the previous sorting cannot be excluded. We do not know how this difference could be assessed. Do you have an idea? But if there was a substantial difference, we would have seen a deviation of the melt around the stake and the immediately surrounding area, which was not the case.
P12 section 4.3 I am a bit confused by the mixed terminology between exposed ice and ice cliffs. At first I thought that with exposed ice, bedsides ice cliffs also patches with very thin or absent debris were meant. However, in this section and also further in the manuscript, it seems like exposed ice and ice cliffs are used interchangeably. Please make this clearer and be consistent with your terminology. If you just mean ice cliffs, I would suggest sticking to that term, since this is most commonly used in literature. We've changed the terminology. We made clear in the methods that the term 'ice cliffs' includes both exposed ice from ice cliffs and exposed ice from water flow channels and then stuck to the term 'ice cliffs' throughout the document.

P12L23
,but the -> ,and This statement has been removed.

P13L1
The promille is a bit confusing and unnecessary here, I think. Just stick to 0.5% and 1.8% in this sentence. New P12 L23 Yes, confusing indeed. Changed as suggested.

P13L10
This also depend on variable rates of ice emergence, which should theoretically be taken into account if an analysis on a subset of a geodetic dataset is performed, also see [Brun et al., 2018]. Ideally, to really look at the effect of ice cliffs and their relative melt, there should be some correction for the downglacier displacement of the ice and cliff. New P12, section 4.5 The objective of this analysis is to understand the potential influence of ice cliffs on elevation change. To do so, we compared area-averaged elevation change with and without ice-cliff and their surrounding areas but could not see substantial differences. Importantly, we do not explicitly compare different time periods. Also, both horizontal displacements and cliff backwasting are inherently accounted for by applying the buffer. Thus, our main conclusions only refer to elevation change per date (not ablation) and are in our view robust and do not require emergence velocities and horizontal displacements. We adjusted the description to make the objective of the analysis and the conclusion drawn from it clear.

P13L26
balances -> stable New P13 L17 Sentence changed to: "… the tongue's elevation barely changed…". P13 section 4.4 Uncertainty ranges should be included in the numbers that are provided in this section. New section 4.6 Done.

P14L7
Pushed down sounds odd. "Travelled downglacier with the flowing ice"? …at a higher rate…, is it relevant if the rate is higher? They are 'pushed' irrespective of the velocity (expect for completely stable ice). We have removed this statement since it doesn't provide any crucial information.

P15 Fig9
Color of the class -0.5-0.5 should be white, in my opinion.

New Fig 10, P14
The problem with white is that it cannot be distinguished anymore from the background, i.e. where there is no data available. We have changed the yellow into yellow-grey to give it a more neutral appearance. Too grey would again move it to the blue category, thus it's a fine line.

P16L6
"at 50 cm resolution" This sub-chapter has been completely removed. The analysis of the difference between the two UAV DEMs was removed because it contained a redundant message. P16L11 remove 'due to higher temperatures' We've removed the sentence, since it was a repetition from above.

P16L14
"There is no clear hint" does not sound very scientific. New P13 L26-28. True. We have strongly changed the text and reduced it to the central content. Also, we incorporated it into the elevation change chapter (4.6), since it is an extension of the same data and also about elevation change patterns.

P16 Fig10
Instead of plotting glacier section number on the x axis it would be more informative to use the mean elevation of a section instead. In that case, you may also consider switching the axes to get an dh/dt gradient kind of plot, similar to figure 17. What does the 'relation' show in the subtable? A regression through sections? Should ideally be filtered for significance. Consider skipping the table entirely. New P15 L1 We agree that when we want to assess the elevation change gradient it makes more sense to look at the changes of elevation-based classes than of spatial sections. Therefore, we have changed the analysis to look at elevation changes of 50 m classes. Additionally, we have reduced the number of periods from eleven to seven. Thus, the additional map as well as the table are not necessary anymore. The whole analysis should now be easier to understand and the figure easier to read.

P16L29
This is not based on your data, right? Provide inline reference. New P9 L5 Probably this is some kind of misunderstanding. Here we present the results of our mass balance analysis.

P18L4
Unclear, please rephrase. New P12 L14 We've rephrased the sentence to: "The displacements from radar interferometry from the 1.5-day period in summer 2017 yield similar results and confirm the quasi-stagnation of the lower tongue …". P18 section 5.1 Could use some subsubheadings New P15 L6 We chose to keep it separated by the different paragraphs. Subheadings would take up additional space and the chapters would mostly be rather short. The chapter has also been shortened considerably (from 913 to 600 words).

P19L9
These uncertainties are a bit unclear to me. A range from about 1 meter to 2 meter. Also, when you convert that to m a-1, doesn't that depend on the variable time span between each observation pair? New P15 L7 Here we provided the range because the uncertainties are different for each DTM. We added a reference to Table 1 where these uncertainties are stated. Yes, the value of the annual uncertainty depends on the length of the time span. We considered this accordingly whenever we provided annual change values.
P19L27-30 I agree. I think it is not fair to compare 1.5-day measurements to those over much longer time spans. There will be quite some intraannual, intraseasonal and probably diurnal variability in velocity. This should clearly be acknowledged prior to the discussion. New P15 L23 We have kept this information at this location where we think it is appropriately placed. Also, the radar information is used as complementary data to get a preciser idea of the velocity pattern and the stagnating glacier part. For this purpose the hint in the discussion to the seasonally restricted meaning is sufficient.

P19L42
Don't miss  here :-) New P15 L30 How could we  We agree that it was an important reference that we missed to cite. P20L8 I don't think the discussed P and T encompass "all glacier-related variables". During the reorganisation process, this statement got removed.

P20 Fig 14
Not a big fan of introducing new figures, data and analyses in the discussion section. Discussion is to discuss the results already presented. We agree. Newly, we included the temperature change information in the mass balance figure. By sharing the same x-axis, the two pieces of information can nicely be linked, even visually.

P20L22
is it a constant increase, or a continuous increase? New P16 L7 We've restructured the sentence to: "At the same time, the debris-covered area increased from 2.8±0.2 km² in 1859 to 3.79±0.01 km² in 1961" P20L28 Exposed-ice? You mean ice cliff? Earlier comment applied to this section as well. We changed the terminology (see earlier comment). It is now 'ice cliff' in this case.

P21L4
Was even moded -> has even moved. Same next sentence. New P16 L21 Done.

P21L7
Attenuated is not the right word here. Attenuate means to make something else smaller, thinner, or weaker. P16 L23 True, thanks for spotting this. We agree that the term is not quite adequate in this case. We've changed it to ".. the retreat rate and area loss gradually slowed down."

P21L14
During…decades -> Over the last 16 years I really don't get these ranges that use the tilde symbol. About 0.80 to 0.98 seems rather specific to me. Why is this an estimate? New P16 L33 Changed to 'In the most recent period (2001 to 2017) …' We tried to keep the discussion on a qualitative level as far as possible. Exact values here would mean the values of three periods including uncertainty ranges, which is long, repetitive, and not really necessary. Therefore we kept the rough values but exchanged the tilde symbol with an 'approx.' P21L15 Same for larger than a range. Isn't >1-2 m the same as >2, essentially? Or do you mean something like greater-than above similar or equal to 2 m? (e.g. ⪎, ≳, ⪆)? New P16 L33 True. We changed it to >1m yr -1 P21L16-17 "As a result length and area changes have been comparably small given the high mass loss". New P16 L35 Thanks for the suggestion! Sentence changed to : "… which might explain the relatively small length and area changes given the high overall mass loss (Figs)." P21L26 "surfacing of debris" During the revision, this statement has been incorporated in other sentences.

P22L17
We always refer to these cavities as voids, as per i.a. [Benn et al., 2012] New P19 L41 We changed 'cavities' to voids as suggested.
P24L26 "A sample" is not very specific. Also why suddenly this new analysis, which does not really provide any novel insights, at the end of the discussion? New P17, chapter 5.4 We have restructured the discussion such that the overview is clear and the sections are linked to respective main messages. This analysis clearly yields that debris-free and debris-covered glaciers can be separated based on elevation change patterns, and that Zmuttgletscher can clearly be grouped to the latter. We have re-structured the discussion to put this analysis in a better position to be used for argumentation. The y-axis of the figure shows only elevation bins and no absolute numbers because the absolute elevation was normalised for comparison. This is necessary to account for the different elevation ranges of the glaciers. To make this clear we will change the axis caption from 'Elevation class' to 'Normalised elevation class'. We also prefer to keep this figure in the discussion since this is the place where we conclude the effect of debris cover. Additionally, this is not a real 'analysis' per se: we have simply used existing DTM data (prepared by Fischer et al. 2015) and plotted it in a comparable way for a selection of glaciers. Thus it would not really fit into the 'results' chapter, and doing so would also lead to additional repetition, which we tried to avoid as much as possible in the revised manuscript. P26L1-… I found the sudden use of procent points a bit strange in this section. Points are a bit irrelevant as they do not show whether something changed from for example no debris to 30% debris, or from 60% to 90%. Surely the actual increases are presented in the original papers? New P17 chapter 5.3 We chose the percentage points because it is easier to compare numbers and the reader doesn't need to calculate, because otherwise it would be necessary to present the debris cover share both at the start and the end of the period. We have now strongly cut down this section, incorporated it into a separate chapter about debris cover changes, and at the same time included the missing information about the absolute debris cover from other studies, which covers a broad range: "These studies found debris cover increases of 2-10 percentage points with a broad range of absolute debris cover from 2-42%." In case you're interested in the absolute numbers, here are the debris cover changes from studies: Mmany of these glaciers show high thinning rates today despite thick debris cover. Due to the lack of observations, most existing studies have neglected the dynamic interactions between debris cover and glacier evolution over longer time periods. The main aim of this study is to reveal this such interactions by reconstructing changes of debris cover, glacier geometry, flow velocities, and surface features of Zmuttgletscher (Switzerland), based on historic maps, satellite images, aerial photographs, and field observations. We show that debris cover extent has increased from ~13% to ~>32% of the total 15 glacier surface since 1859, with decadal and that today the debris is sufficiently thick to reduce ablation compared to bare ice over much of the ablation area. Despite the debris cover, the glacier-wide mass balance volume loss of Zmuttgletscher is comparable to that of debris-free glaciers located in similar settings, whereas changes in length and area have been small and delayed in comparison. Increased ice mass input in the 1970s and 1980s resulted in a temporary velocity increase, which led to a local decrease in debris cover extent, a lowering of the upper margin boundary of the ice-cliff zone, and a strong 20 reduction of in ice-cliff area, indicating a dynamic link between flow velocities, debris cover, and surface morphologys.
Since ~20015, the lowermostst ~1.5 km of the glacier have beenare quasi-stagnant, despite a slight increase in surface slope of the glacier tongue. We conclude that the long-term glacier-wide mass balance is mainly governed by climate. tThe debris cover governs the spatial pattern of volume loss elevation change without changing its glacier-wide magnitude, which is due to we explain by the large extended ablation area and strong the enhanced thinning in regions with thin debris further up-25 glacier and in the regions of areas with abundant meltwater channels and ice cliffs. At the same time rising temperatures lead to increasing debris cover and decreasing glacier dynamics, thereby slowing downattenuating length and area losses.

Motivation and Oobjectives
Debris-covered glaciers have been observed to show a delayed adjustment of their length to climatic changes (e.g. Ogilvie, 1904;Scherler et al., 2011;Banerjee and Shankar, 2013). This behaviour can be explained by melt reduction due to 30 insulation by the debris layer, which commonly increases in thickness towards the terminus (Östrem, 1959;Nakawo et al., 1986;Nicholson and Benn, 2006;Anderson and Anderson, 2018), and is expected to distinctly prolong the glacier's response time (Jóhannesson et al., 1989). In some regionsSeveral studies showed that debris-covered glaciers exhibit similar thinning rates of volume changes asto debris-free glaciers in the Himalayas (e.g. Kääb et al., 2012;Gardelle et al., 2012;Nuimura et al., 2012). Explanations pProposed reasons for this behaviour are:include the emergence of surface features -35 typical for debris-covered glaciers -with exposed ice (e.g. ice cliffs, water flow channels, and ponds), stronger enhanced thinning further up-glacier that compensates for the debris-induced melt reduction, reduced mass flux from the accumulation areas and decreasing emergence velocities at the tongue, and englacial ablation with subsequent roof collapsing of the glacier surface (e.g. Pellicciotti et al., 2015;Vincent et al., 2016 ;(Vincent et al., 2016); Ragettli et al., 2016;Banerjee, 2017;Brun et al., 2018). (Brun et al., 2018)). Ice cliffs and ponds can enhance ablation in comparison to debris-covered surfaces and even debris-free ice and are common features on debris-covered glacier tongues (Benn et al., 2012;Brun et al., 2016). Especially during periods of negative mass balance, the down-glacier increase of in debris -cover thickness reduces ablation through its insulating effect and can lead to a lower -and even reversed -mass balance gradient (Nakawo et al., 1999;Benn and 5 Lehmkuhl, 2000;Benn et al., 2012;Ragettli et al., 2016;Rounce et al., 2018). Over time, this reduction in ablation can lead to a decrease in surface slope and, consequently, driving stress andOver time, the surface slope of the glacier tongue therefore often reduces, which can lead to a decrease in driving stress and thus ice flow velocity (Kääb, 2005;Bolch et al., 2008b;Jouvet et al., 2011;Rowan et al., 2015;Dehecq et al., 2019). Furthermore, with an increase in equilibrium line altitude (ELA) the englacial debris melts out fasterearlier, leading to an extended debris cover further up-10 glacier (Benn et al., 2012;Kirkbride and Deline, 2013;Carturan et al., 2013). As a resultconsequence of reduced ablation and driving stress, strongly heavily debris-covered glaciers often have long, flatgently-sloping, and low-lying tongues with low flow velocities or even stagnant parts at their tongues (Benn et al., 2012). Current research on many debris-covered glaciers is mostly focussed on investigates the involved processes, i.e.such as ablation under abeneath the debris cover and in areas of ice cliffs and ponds, ablation and thinning at of the glacier tongues, surface flow velocities, or changes of in debris 15 cover over time, or surface flow velocities, at many debris-covered glaciers, specifically in the Himalaya (e.g. Hambrey et al., 2008;Bolch et al., 2012;Benn et al., 2012;Dobhal et al., 2013;Ragettli et al., 2016;Pellicciotti et al., 2015;Gibson et al., 2017).; A few some studies have started to use integrate the existing understanding and interactions knowledge about the mutual influence of some of these processes into numerically models for the evolution of such glaciers' development (Jouvet et al., 2011;Rowan et al., 2015;Anderson and Anderson, 2016;Banerjee, 2017). 20 However, most studies face difficulties leading to persistent shortcomings in our understanding of the development of debris-covered glaciers: (i) t: time series are often too short to inform on the full and with a few decades well below the expected response time of larger glaciers;; (ii) investigations are often local because repeated tongue-wide data are sparse, especially for longer periods;; (iii) by considering only one or few variables are considered, and thus the mutual influences of e.g. changing debris cover and glacier geometry cannot be assessed; ; (iv) studies at glacier-scale over more than a decade 25 have mostly been conducted in the Himalaya (e.g. Bolch et al., 2008b;Ragettli et al., 2016;Lamsal et al., 2017).
To better understand how a changing debris cover affects glacier geometry, flow velocities, and surface features, and how it is in turn affected by these variables, it is necessary to consider the long-term development of glaciers beyond their potential response times. Few studies have investigated the temporal evolution of debris cover on glaciers (overview in Kirkbride and Deline, 2013), or the evolution of debris-covered glaciers over time (e.g. Agata & Zanutta 2007;Capt et al., 2016). . Several 30 studies observed an increase inof debris cover on glaciers during negative mass balance periods (e.g. Bolch et al., 2008a;Quincey and Glasser, 2009;Kirkbride and Deline, 2013).; Quincey and Glasser, 2009;Bolch et al., 2008a). Because of the overall negative mass balance trend, glaciers are and will be increasingly affected by debris cover. It is important to understand how the magnitude of this increase and how it affects the geometry and dynamics of the glaciers in order to simulate their future development (Jouvet et al., 2011;Anderson and Anderson, 2016). 35 Learning about the long-term effects of debris cover will be valid not only for glaciers in the Alps, but for glaciers worldwide that undergo general retreat. Large debris-covered glaciers in the Himalaya and Karakoram -which are in the centre of attention due to their importance -are not ideal candidates for long-term investigations due to their long response times (>50 years) and data scarcity before the 1960s. In contrast, the long history of length monitoring as well as the availability of topographic maps and aerial photos from the mid-19 th century onwards make glaciers in the Alps and 40 especially in Switzerland Swiss glaciers well suited sites to study long-term developments;. tIn Switzerland, the earliest maps (1850s, 1870s) already include debris symbols and contour lines, and stereo imagery throughout the 20 th century allows for detailed 3D surface investigations.
In this study we aim to understand the long-term geometric evolution and dynamics of debris-covered glaciers at the examplethrough the study of Zmuttgletscher in the Swiss Alps. This medium-sized valley glacier has been going through the transition from a mostly debris-free glacier in the late 1850s to one that is almost completely debris-covered in its ablation area today. WWe quantify the evolution of geometry (length, area, elevation, slope), mass balance, flow, and debris cover at a high spatial and temporal resolution since the end of the Little Ice Age (LIA) around 1850. We use these data to investigate 5 the relation and interactions between geometry evolution, ice flow, and debris covecover and the related driversr. We further analyse the occurrence of ice cliffs and their role for the long-term glacier evolution. Zmuttgletscher is located in a relatively dry region at the main divide of the Alps, thus it receives precipitation from both northern and southern weather situations. There are no direct measurements at higher elevations, but model estimates suggest 15 values between 800-1500 mm (MeteoSwiss, 2014). Glaciological mass balance measurements at the nearby debris-free Findelgletscher (15 km distant) suggest end-of-winter accumulation around 0.8-1.5 m water equivalent (w.e.) (Sold et al., 2016). However, at Zmuttgletscher, avalanching additionally contributes to accumulation. When including contributing rock walls and lateral moraines, the total area available for accumulation is up to 22 km² (restricted to areas >30° slope angle,  (Figure 1). At present, the glacier has a surface area of ~16 km² with a substantial part of debris cover in its ablation area, originating from the surrounding rock walls.
Zmuttgletscher lies in the area of the Dent-Blanche nappe, which is part of the East-Alpine nappe. Most summits and rock walls consist of slightly metamorphosed crystalline and magmatic rocks (various types of gneiss, gabbros and diorites (Gouffon and Bucher, 2003). 40 4 Zmuttgletscher is located just north of the divide of the main Alpine ridge. The Rhone and Aosta valleys a few kilometres to the north and south, respectively, are well shielded by high mountains to the north and west and belong therefore to the driest regions in the Alps with precipitation values around 600 mm (Grächen, Mattertal: 653 mm;MeteoSwiss, 2017). Precipitation is mostly linked to northern and southern weather situations with a strong western component. No direct measurements at higher elevations are available; model estimates suggest values between 800-1500 mm in the area of Zmuttgletscher 5 (MeteoSwiss, 2014). Glaciological mass balance measurements at near-by almost debris-free Findelgletscher show end-ofwinter accumulation around 0.8-1.5 m water equivalent (w.e.) (Sold et al., 2016) and similar values can be assumed for Zmuttgletscher. However, avalanching importantly contributes to accumulation. Including contributing rock walls and steep lateral moraines, the total area available for accumulation that is eventually being transported onto the ice is up to 22 km² (restricted to areas of >30° slope angle). The glacier was part of the regular length monitoring programme of the Swiss glacier monitoring network (GLAMOS) until 1997, when it was abandoned due to the uncertain position of the debris-covered, stagnant terminus. Furthermore, dendrochronological analysis tried to reconstruct former historic advances (Schneebeli and Röthlisberger, 1976), and an early study by Haefeli (1952) investigated ice deformation.

Data and Mmethods 5
The analysis of this study was based on input fromOur analysis are based on several sources: topographical maps, airplane-based and UAV-based aerial images, and (Table 1), satellite images (Table 1), in addition to various field observations, and long-term air temperature measurements. The use of these data is discussed in the respective sections below. Areal images were available from (i) post-war mapping flights by the American military, (ii) national flight campaigns for topographic map production, (iii) specific flight campaigns for glacier monitoring purposes conducted by 10 Swisstopo, (iv) fixed-wing UAV flights (using a SenseFly eBee Classic; SenseFly SA) in 2016 and 2017, and (v) Pléiades tri-stereo images from 2017.

Generation of digital terrain models (DTMs) and orthophotos generation
Glacier surface information was extracted from stereo-photogrammetric DTMs generated from aerial images and the Siegfried map (Table 1) Pléiades tri-stereo images from 09/2017 • To obtain a DTM representing the glacier tongue in 1879 the Siegfried map from 1880 (Swisstopo, 2018a) was 5 georeferenced using ground control points (GCPs). Subsequently, the contour lines were semi-automatically digitized by separating the differently coloured (blue = on-glacier, brown = off-glacier) contour lines from other symbols (Siedler, 2011).
Their elevation information was extracted and interpolated using ArcGIS´ Topo-to-Raster tool. The point of origin for elevation measurements, situated at the shore of Lake Geneva, changed in the 1930s from 376.2 m to 373.6 m (Swisstopo, 2018b) and the elevations derived from the Siegfried map were corrected accordingly. 10 The The time series of DTMs from 1946 to 2017 (apart from 2010 and the Pléiades DTM from 2017) wasere created with photogrammetric methods using Structure-from-Motion software packages (Agisoft LLC, 2016;Pix4D, 2016).; cf. Mölg and Bolch, 2017). For each dateDuring this process, a point cloud was produced from the available number of input images (4-29) and was then georeferenced using a set of 10-15 GCPs, i.e. reference points in the images that could be referred to points on stable ground, taken from Swissimage 2013. The quality of the DTMs is comparable to DTMs from traditional 15 photogrammetric software (cf. Mölg and Bolch, 2017). The quality of the DTMs is comparable to DTMs from traditional photogrammetric software (Mölg and Bolch, 2017). Their uncertainty in elevation -defined byas the standard deviation over stable terrain around the glacier tongue -mostly lies within 2 m depending on the resolution of the aerial images, their number and image quality factors (somewhat higher values were obtained for the DTMs from 1879 and 1946; Table 1). The uncertainties were derived for terrain with comparable steepness to the glacier surface (<30°) and are assumed constant in 20 space, although DTM quality decreases in steep areas (e.g. rock walls surrounding the glacier). The DTM from 2010 (SwissAlti3D) was taken as a final product from Swisstopo and had also been produced from aerial stereo images; it has a nominal resolution of 2 m and a vertical accuracy of 2 m (Swisstopo 2018). In addition, we generated a glacier-wide DTM for 2017 (spatial resolution 1 m) using Pléiades tri-stereo, high-resolution satellite imagery, from Centre national d'études spatiales -Pléiades, with a ground resolution of 1 m. The tri-stereo image was produced based on photogrammetric 25 principles using PCI Geomatica (Geomatica 2016). All DTMs were co-registered before further analysis by using the analytical approach by Nuth and Kääb (2011), followed by a second-order trend surface correction to eliminate remaining elevation differences (Pieczonka et al., 2013).
Orthophotos were generated by rectification of the stereo images using the respective DTM. These images were subsequently used for delineating the glacier boundary and mapping of exposed ice. We further used the Swissimage by 30 Swisstopo -orthorectified image data from aerial imagery -that was available for the years , 2007(Swisstopo 2018.

Glacier area and length
Glacier areas waswere extracted measured from the Dufour map (1859), the Siegfried map (1879), and all available orthophotos from this study and from Swisstopo (2005Swisstopo ( , 2007Swisstopo ( , 2010Swisstopo ( , 2013Swisstopo, 2010; Table 1) by manual  35 digitizationing. The Dufour map (map sheet 22, section 8, number 485) from 1859 is was the first map containing elevation information (in the form of contour lines) and distances acquired with modern methods (Wolf, 1879;Graf, 1896;Rastner et al., 2016). The extent of the glacier and the supraglacial debris could be extracted from the map, whereas its elevation information was disregarded due to strong, non-linear, horizontal distortions. The Siegfried map in 1880 was a follow-up, containing elevation information at Zmuttgletscher from 1879 (Swisstopo, 2018a). The map was georeferenced using ground 40 control points (GCPs), i.e. reference points indicated in the map that could be referred to stable points today. The hanging glaciers at the north face of Dent d'Hérens have been included in the glacier area, since they contribute mass to the main glacier through roughly regular 'ice fall' events. This time series of maps and orthophotos resulted in a glacier area information value for each corresponding date since 1859. The mapping quality is commonly often lower in debris-covered compared to debris-free areas (Paul et al., 2013), but the high resolution of the images allowed correct interpretation ofs for correctly interpreting the glacier margin. The glacier boundary in the accumulation area of STG to Glacier de Ferpècle and 5 Haut Glacier de Tsa de Tsanto the West was taken from the 2010 ice divide and was kept stable constant over time. The hanging glaciers at the north face of Dent d'Hérens have been included in the glacier area, since they contribute mass to the main glacier through frequent ice avalanches.
Front variation measurements were conducted by using the glacier outlines for each date. Along the ice front -perpendicular to the flow -the change was measured at distances of 100 m and then averaged (Koblet et al., 2010;Bhambri et al., 2012). 10 For the comparison to other glaciers we used GLAMOS length variations (GLAMOS, 2018), which were acquired in the field with the same concept of several parallel measurements, equidistant by 50 m. At Zmuttgletscher, GLAMOS measurements have regularly been conducted almost annually from between 1892 until and 1997, but with serious data issues before 1946.an additional data point was added for 2010. TThe GLAMOS data for Zmuttgletscher since 1946 are almost identical to our own relative length change record ( Figure S1) and are therefore not further used. To properly interpret 15 the observed length changes information and the potential influence by debris cover, we compared the variations from our measurements to the ones by GLAMOS, data of and ten other finally compared length variations over time from various Swiss glaciers.
The uncertainty of both area and length results are estimated to lie within ± 1 / 2 pixel (Table 1) along the glacier boundary or at the start and end of the length profile (Hall et al., 2003;Granshaw and G. Fountain, 2006;Bolch et al., 2010;Bajracharya 20 et al., 2015). For the calculation of the changes, the uncertainties of the two respective dates are combined by error propagation, analogue to the mass balance data (Hall et al., 2003;Zemp et al., 2013).

Debris cover, on-site ablation and on-site air temperature
Debris cover extent for the orthophotos and historic maps was manually digitized using the orthophotos as well as the two historic maps which already contain a debris cover symbol ( Figure 2a, and b). Further, the debris extent of the Siegfried map 25 (1879) could be verified using two photographs taken in 1894 (Figure 2c). This information is valuable to limit the debris extent at and up-glacier of the confluence of TMG and SBG, which is the region of the strongest changes. A continuous debris cover has reflection characteristics different from ice and can be distinguished visually. We mapped only areas of continuous debris cover, thus avoiding the inclusion of sparse debris coverThe transition from debris-free to debris-covered is continuous, and regions of sparse debris cover -which might even lead to increased ablation -were not included. We 30 mapped only coherent areas of continuous debris cover, thus avoiding to include mixed pixels that indicate a sparse debris cover. The debris extent of the Siegfried map (1879) was verified using two photographs taken in 1894 (Figure 2c). This information was valuable to limit the debris extent at and up-glacier of the confluence of TMG and SBG, which was the region of the strongest changes ( Figure 2d). Because of varying properties of the orthophotos (RGB vs. panchromatic, resolution, contrast, texture) an automated method would not deliver equally satisfying results and was not applied.  (Fahnestock et al., 1992;Scambos et al., 1992;Conrad et al., 2015). The results were filtered using a visually defined threshold of correlation quality. Outliers were manually removed, e.g. in the area of strong ice-cliff change or pro-glacial water surfaces.
Between 22 nd and 24 th of August 2017 the tongue of Zmuttgletscher was observed with a terrestrial radar interferometer 5 (GPRI) developed by GAMMA (Werner et al., 2008). The GPRI was installed on a hill about 3 km away from today´s terminus (Figure 1). Measurements were acquired every minute for 1.5 days with a final range resolution of 3.75 m and an azimuth resolution of 7 m at a slant range of 1 km. The interferograms were determined with a standard workflow following Caduff et al. (2015) using the Gamma software, were stacked over a window of 8 hours to reduce noise, and were afterwards unwrapped using a stationary point on the ground. The unwrapped phases were then converted to line-of-sight displacements 10 according to Werner et al. (2008), whereby negative velocities were considered to be noise and filtered out. The results were georeferenced by rotating the displacement map to best match with the DEM25 (Swisstopo, 2005). Afterwards, the data was resampled into the new grid using nearest neighbour interpolation. To assess the uncertainties in the velocities we looked at the difference from zero in measured displacement of 10 stationary points. This results in an uncertainty of the stacked velocity maps of ±0.03 m/day. 15

Surface features
Ice cliffs, exposed ice at supraglacial meltwater flow channels from supraglacial melt streams with exposed ice, and supraglacial lakes w were extracted in a semi-automatic process using an object-based approach with Trimble eCognition  (Sakai et al., 2002;Buri et al., 2016). These changing slope areas were often separated by segmentation, and were in a second step manually selected and combined into one polygon per ice cliff or lake. Supraglacial channels cover only small areas and are thus incorporated into the term 'ice cliffs'.Using the approach above, water flow 25 channels and ice cliffs could not be mapped as separate features and are summarised under the category 'exposed ice area'.
The described approach is effort efficient and assures a low uncertainty level, which is in the order of ±½ pixel.
In order to assess the importance potential influence of exposed-ice areascliffs with respect to on surface elevation change, we investigated their changes in area at the lower section of on the glacier tongue over timefor different time periods. A 35 m buffer was generated around the exposed-ice-cliff polygon for both date 1 and date 2, and the overlapping area of these 30 two buffer zones defined the category of exposed-ice-cliff areas. For each period the average surface elevation change in this category was compared to the surface elevation change at the rest of the tongue (yellow area in Figure 9).

Reconstructing the glacier surface topography for the 19 th century
At the end of the LIA many Alpine glaciers reached a maximum Holocene extent, similar to the LIA maximum around 1650, and accordingly moraines were built up or former high stands were again reached, approximately in the 1850s (Steiner et al., 35 2008;Zumbühl et al., 2008). Moraine elevation is available for large parts around today's ablation area of Zmuttgletscher (indicated in Figure 1). The maximum elevation of lateral moraines was used to interpolate the glacier surface elevation using the Natural Neighbour method (Sibson, 1981); we assumed a plane across the glacier surface due to the small curvature indicated in the Dufour map (Graf, 1896;Swisstopo, 2018b). This LIA surface elevation was not used for calculating geodetic mass balance because of the incomplete areal coverage. 40

Formatted: Not Superscript/ Subscript
The contour lines of the Siegfried map were semi-automatically digitized by separating the differently coloured (blue = onglacier, brown = off-glacier) contour lines from other symbols (Siedler, 2011). Their elevation information was extracted and interpolated using ArcGIS´ Topo-to-Raster tool to achieve a distributed DTM across the glacier tongue. The point of origin for elevation measurements, situated at the shore of Lake Geneva, changed in the 1930s from 376.2 m to 373.6 m (Swisstopo, 2018b) and the elevations derived from the Siegfried map were corrected accordingly. Elevation information 5 from the Dufour map was disregarded due to strong, non-linear, horizontal distortions.

Surface elevation changes and geodetic mass balance
The time series of surface elevation change data over the glacier tongue was restricted to the overlapping area of all available DTMs (except the moraine-derived LIA DTM), from 1879 to 2017 (11 dates), and which reaches up to ~2750 m a.s.l. and covers >50% of today's ablation area of TMG as well as the lowest parts of STG and SBG~20% of today´s total glacier area. 10 The surface elevation change maps were generated and investigated for each individual period as well as for the entire time Glacier-wide Ggeodetic mass balance estimates wasere calculated for periods between dates for which where the DTMs covered most large parts of the glacier surface (1879,1946,1977,1988,2001,2010,2017). In areas of data voids or artefacts 20 in the DTM, especially in higher elevations, Due to DTM insufficiencies (artefacts, voids), certain areas were not covered, especially in higher elevations, and no surface elevation change values could be calculated there. In order to reduce the sensitivity to data voids and outliers, 100 m elevation bins (starting from 2100 m) were used for the elevation change calculation. In elevation bins without DTM coverage, the average of three extrapolation methods was used as the final mass balance value. The following three methods were used to fill the data voidsThese missing values were filled using a range of 25 methods: (1) a linear relationship between elevation and elevation change (e.g. Kohler et al., 2007;Kääb, 2008);, which is based on the strong respective correlation up to R²=0.93); (2) using the same relationship but additionally setting the elevation change values to 0 when they become positive in the highest elevations; (3) using the mean elevation change value of the uppermost elevations that contain data and applying it to the glacierised area above. Independent of the method The area-weighted average elevation change of all bins was it was necessary to work with elevation bins of 100 m (starting from 30 2100 m) to get a representative value on which to base the calculations. Reason for this was the incomplete coverage of some glacier parts and thus a certain susceptibility to outliers, which could be minimised by the use of elevation bins. The average elevation change per elevation bin was multiplied with the area share of this bin of the total glacier area and summed up to reachto calculate the average elevation change of the total glacier surface. We assumed an average density of 850 ±60 kg/m³ ((Sapiano et al., 1998;Huss, 2013) to convert the volume to mass changes. 35

Uncertainties
The total uncertainty of the surface elevation estimates and the mass balance consists of in (i) the accuracy of the individual DTMs, (ii) the filling of the empty elevation zones, (iii) the glacier volume density changes and density conversion, (iv) the debris volume changes, and (v) the DTM co-registration. Ice density changes are on average assumed to be negligible over a longer time span (Huss, 2013) but we . As no information is available about changes in snow and firn volume, we assumed 40 zero change for that, but included an uncertainty of ±60 kg/m³ for the density to mass conversion (σ dens ), cf.analogue to Huss Formatted: Not Superscript/ Subscript (2013). Debris volume changes were not accounted for because of lack of such data, but they would have a negligible effect on the total uncertainty in volume change. Because not all DTMs used for mass balance calculation covered the entire glacier surface, the missing surface area had to be filled with values from another point in time, in this case the reference DTM from 2010. This procedure introduces an uncertainty since the area of the elevation bin might have changed over time, but the order of change is small due to the absence of big topographic steps in the glacier hypsography (Figure 1c). The error from 5 the filling of the elevation zones is taken into account as the standard deviation between the total mass balance after applying the three different interpolation methods (σ fill ). This uncertainty measure was combined with the uncertainties of the two DTMs framing the respective period (σ DTM1 and σ DTM2 , standard deviation over stable terrain) to be used as the total uncertainty for the resulting mass balance by applying the law of error propagation (e.g. Hall et al., 2003;Zemp et al., 2013) and hence given by 10 Uncertainties range from ±0.12 to ±0.2 m w.e. /yr -1 , whereas the uncertainty for the total period 1879-2017 is ±0.03 m/yr yr -1 (Table 2). .

Surface elevation change for different debris thickness classes 15
On the time scale of several years, the effect of various debris thicknesses may impact on surface elevation changes. Using the debris thickness observations from field measurements, we defined debris thickness classes of thin debris (~5 cm) and thick debris (>15 cm) as well as of bare ice. We mapped representative areas of these classes for several periods between 1995 and 2017 (1995-2001, 2001-2005, 2005-2010, 2016-2017) and within these areas of ~20.000 m² we calculated the average elevation change values for each period. 20

Surface flow velocities
Surface flow velocities provide important information about the glaciers' dynamical state and its change over time.
Automatic feature tracking methods were not feasible for the time before 2005 because the time differences and, thus, displacements of the features were too big. Therefore we manually tracked boulders to infer flow velocities along the debriscovered part of the tongue of Zmuttgletscher as well as on the lower branch of SBG. The tongue was divided into eleven 25 sections according to differences in dynamic state and ice flow units. The individual measurements were averaged for every time period and section, respectively, to achieve a comprehensive picture of the dynamic changes ( Figure 12). For the periods 2005-2007 and 2016-2017 we extracted flow velocity fields using the IMCORR module in SAGA GIS (Fahnestock et al., 1992;Scambos et al., 1992;Conrad et al., 2015). The results were filtered using a visually defined threshold of correlation quality ('Strength'). Outliers had to be manually removed, e.g. in the area of strong ice cliff change or pro-glacial 30 water surfaces.
On 22 nd -24 th of August 2017 the tongue of Zmuttgletscher was observed additionally with a terrestrial radar interferometer (GPRI) developed by GAMMA Remote Sensing and Consulting AG, Switzerland. The GPRI was installed on a hill about 3 km away from today´s terminus (Figure 1). Measurements were acquired every minute for 1.5 days with a final range resolution of 3.75 m and an azimuth resolution of 7 m at a slant range of 1 km. The interferograms were determined with a 35 standard workflow following Caduff et al. (2015) using the Gamma software. The interferograms are stacked over a window of 8 hours to reduce noise and afterwards unwrapped using a stationary point on the main land. The unwrapped phases were then converted to line-of-sight displacement according to Werner et al. (2008), whereby the negative velocities are considered as noise and filtered out. The results were georeferenced by rotating through an estimated angle based on the best match with the DEM25 (Swisstopo, 2005). Afterwards, the data was resampled into the new grid using nearest neighbour 40

Formatted: Superscript
Formatted: Superscript 13 interpolation. To access the uncertainties in the velocities we looked at the difference from zero in measured displacement of 10 stationary points. This results in an uncertainty of the stacked velocity maps of ≤ 0.03 m/day.

Geodetic mass balance
Zmuttgletscher's long-term mass balance from 1879 to 2017 has been negative at -0.31±0.03 m w.e. yr -1 ( Table 2). The mass 5 balance was least negative between 1879 and 1946 (Figure 3a). Also during the climatically favourable period in the 1970s and 1980s Zmuttgletscher yielded a relatively small mass loss, even below the Swiss-wide mean. Between 1977 and1988 the glacier-wide mass balance stayed negative (-0.27±0.18 m w.e. yr -1 ), although certain areas on the tongue showed an increase in elevation (Figure 10c, d, e). After 1988, the mass balance became more negative again, even though the lowest areas on the tongue still showed stable or even increasing surface elevation (Figure 10f). The negative trend continuously 10 intensified and after 2001 Zmuttgletscher's mass balance has become more negative than the Swiss-wide mean ( Table 2).      Generally, the debris cover extent has expanded up-glacier in TMG and SBGall parts of the glacier. The extension was very pronounced along the surface of SBG and in the central part of the main glacier tongue (Figure 6). In both areas the debris 5 extent has expanded to above the ice fall into the former accumulation areas and starts now close to the foot of the contributing rock walls of Dent d'Hérens (TMG) and Dent Blanche (SBG). The debris extent also cover grew stronglystrongly expanded also at the glacier margins below Stockji, due to further input from moraines and the disconnection of a contributing tributaryies. By 2010, the ice fall at STG had thinned sufficiently to disconnect and expose the rock wall beneath, from which a small rock fall detached between 2010 and 2013 ( Figure 6)Below the ice fall of STG, a 10 small rock fall occurred between 2010 and 2013 (indicated in Figure 5) after the exposure of the rock wall, covering an area of approx. 170x300 m. This debris mound covered an area of approx. 170x300 m and is now being slowly transported downglacier and , while spreading laterally due to preferential ablation and fluvial transport.

Debris cover characteristics
Field investigations revealed a homogeneous debris cover in some regions with superficial stone sizes in the top-most layer mainly between 10 and 30 cm in diameter, and a much more heterogeneous debris cover in other regions, with pebble-sized stones to metre-scale boulders ( Figure S4) in the metre scale (). In most regions the debris consists of more than one 'layer' (i.e. stones lie on top of each other). The typical base layer of the debris In such cases consists of fine-grained material (sand 10 and even silt), and is overlain by a few centimetres of pebbles.the debris is sorted: fine-grained material (sand and even silt) lies directly on the ice, followed by a few centimetres of pebbles; above these typical base layers follow larger stones without any more specific sorting. The rock walls' different geology leads to differently coloured debris material, which has been deposited on the same location and led to elongated, coloured debris bands along the glacier tongue. The thickness of the debris layer along the each transect varieds from below 5 cm to over 70 cm. The thickest areas were found on the 15 elongated ridge, especially on the steep, southern slope (Figure 7a, transect length 400-500 m). The average thickness along the upper transect was 16.3 cm ±1.3 cm.

10
Ablation measurements from seven locations and over several periods from August 2016 to October 2018a seven weeks period in summer 2017 show an 'Östrem-like' behaviour with respect to debris thickness (Figure 7b). These results, from stakes at elevations between 2300 m and 2600 mdata indicate that melt is strongly dependent on debris thickness and much less on elevation. Compared to the reference stake at 2606 m on debris-free ice, field measurements in summer 2017 (7 weeks) showed a reduction of ~15% and ~50% for a debris thicknesses of 10 and 25 cm, respectively. Our results yield 15 aThese values are somewhat smaller melt reduction by debris compared to than in the literature values, where melt reductions of ~40-60% and 70-80% for ~10 cm and ~25 cm, respectively, were found (Östrem, 1959;Mattson et al., 1993;Nicholson and Benn, 2006;Brock et al., 2007). Note that our debris-free reference stake is at the highest elevation of all stakes, which may partly explain our underestimation in melt reduction by debris. The backwasting rate of ice cliffs was 1.3 times than at the debris-free reference stake. The uncertainty of the ablation measurements is ~3 cm, which results in an ablation rate uncertainty in the order of ±0.09 cm d -1 .

Surface flow velocities
Overall, there is a trend of decreasing velocities since the first measurement period 1961-1977 until today, but with a clear 5 phase of acceleration in the late 1970s and 1980s. Velocities in the lowest section of the glacier tongue (yellow) as well as the Lower SBG have decreased from 10-20 m yr -1 in the 1980s to almost zero (<3 m yr -1 since 2005) and can be considered close to stagnant. The central section (green) showed values of 30-40 m yr -1 until the 1980s; afterwards velocities strongly decreased to ~5-10 m yr -1 today. Velocities in the highest section (blue) have decreased over the same period from 50-60 m yr -1 to ~15 m yr -1 . In the periods 1977-1983 and 1983-1988 a velocity increase of close to 50% was observed in all areas, 10 with a slightly delayed signal down-glacier. It was followed by a strong velocity decrease in the 1990s that became weaker thereafter. In all time periods central glacier ('mid-left', 'mid-right') parts moved faster than the margins, and velocities strongly decreased down-glacier (Figure 8; Figure S5). Uncertainties are in the order of ±0.1 m yr -1 .

Ice cliffs and related elevation change
Even though the glacier has become more and more debris-covered, we did not find any clear long-term trend in area and location of ice cliffs For the location and area of exposed ice we did not find any clear trend in number, area or location 25 (Table 3). Areas of exposed ice They are mostly almost exclusively located in the lowerst parts of the glacier tongue until the terminus (yellow area in Figure 9). The area covered by exposed iceThe total ice-cliff area decreased by over 40% in the 1970s and early 1980s but but increased again afterwards. The average slope of the exposed ice areas extracted from the DTMs is ~26° with no trend over time, but the variations are rather a result of variable DTM resolution and quality. In contrast, own field measurements show slopes of 40-50° for the ice cliffs. Due to this relatively low slope (in comparison to 30 field measurements and literature, e.g. Reid and Brock, 2014) and the small total area, also the tThe topographically corrected exposed-iceice-cliff area amounts to less than 0.5%5‰ of the total glacier area, and less than 1.8% of the debriscovered glacier area (Table 3). As only Only very ffew small ponds have been detected (<5 per time step) and due to their small number and area they have not been further analysed. Over the total entire period, the average elevation change over the exposed-iceof ice-cliff areas (cliffs + buffer overlap) was 1.2 times higher than over the rest of the tongue (Figure 9). During the 1970s and 1980s the thinning rates were, however, almost the same and close to zero whereas since 1988 the ratio of average elevation change ranged the factor of increased elevation change around exposed-ice areas was between 1.5 and 1.7 consistently higher than on average. When considering 10 average thinning rates over the entire glacier tongue after 1995, the inclusion of exposed-ice ice-cliff areas enhances the surface lowering only by less than 5 %, and thus the presence of ice cliffs does not seem to affect the mass loss substantially.

Tongue-wide surface elevation changes
Since the end of the LIA, the surface elevation of the Zmuttgletscher tongue below 2750 m a.s.l. has almost continuously lowered ( Figure 10). Between 1879 and 1946, the average elevation change was negative at -0.63 ±0.13m yr -1 , and substantially lowered to -1.75 ±0.34m yr -1 between 1946 and 1961. In the 1970s and 1980s the tongue partially thickened, 5 which is most pronounced at the tip of STG tongue (+1.2±0.29 m yr -1 from [1977][1978][1979][1980][1981][1982][1983] and in the main trunk of the tongue, while the terminus area was still substantially thinning. Between 1977 and 1983 the tongue's average surface elevation barely changed, at a change rates of -0.11 ±0.29m yr -1 . Between 1983 and 1995 the upper part of the tongue already started to thin again whereas in the lowest 1 km a mixed signal of thickening and thinning occurred (Figure 10f), consistent with a slight advance at the centre of the terminus after 1995. After 1988 the tongue began to thin substantially at a rate of elevation 10 change turned with -1 m yr -1 more negative, becoming even more negative with <-2 m yr -1 after 2005 (see Figure S16  High local values of elevation decrease are observed in regions with high proportions of exposed ice or where no debris cover was present (panels f, h, I in Figure 9). Panels d and e of   Figure S7 to S16 for enlarged panels. .

Surface elevation change for different debris thickness classes
The average total elevation change for the period 1995-2017 was -13.6 m for bare ice, -6.8 m for thin debris, and -4.5 m for thick debris (uncertainty of ±1.1 m), which corresponds to thinning rates of 1.1 m/yr, 0.57 m/yr, and 0.38 m/yr, respectively.
The relation of thick debris to bare ice ranged between 0.29 and 0.63 (average = 0.33). It is important to note that the different debris thickness classes are -by nature -located in different elevations and, hence, do not share the exact same 5 climate.
With the UAV-derived DTMs from September 2016 and October 2017 the surface elevation change of one year in a resolution of 50 cm was calculated for the same glacier extent as displayed in Figure 9. The average elevation change in the overlapping area over the 13 months is -2.4 ±1.8 m but is very heterogeneous. It reveals along-flow structures in elevation change as well as spatial variation patterns related to ice cliffs and flow channels. 10

Evolution of the mass balance gradient
The pattern of the elevation change gradient has strongly shifted over time . A glacier mass balance gradient is typically inclined towards lower elevations due to higher temperatures, whereas the insulation effect of a continuous debris layer can flatten this gradient. To estimate the strength of this debris cover effect, we divided the main tongue of Zmuttgletscher into several sections from 1 (terminus position of 1946) to 24 (highest glacier elevations of the overlapping DTMs, which are still 15 to some extent debris-free today; Figure 10b). There is no clear hint of an increasing balancing effect of debris cover over time (Figure 10a). A reversed gradient existed in the period 1983-1995, due to the 'wave' of additional mass flux travelling down-glacier, which is also the reason for the less negative and partly even positive values.

b): glacier sections from 1 at the terminus of 1946 to 24 at ~2700 m (today). c): Average elevation change values (in metres per year) of the five highest and lowest sections per period.
The average elevation change rate of the lowest and the highest five sections (disregarding the lowest sections after melt-out) yield the strongest contrasts; the relation between these two regions has changed over time and has been above 1 in recent years (Figure 10c), even though the elevation change has reached high absolute values. We can discern a levelling of the 5 elevation change in the lower part of the glacier tongue, which we attribute to the effect of debris cover, although the two lowest sections show still the strongest thinning.

Glacier mass balance
Zmuttgletscher's long-term mass balance from 1879 to 2017 has been negative (-0.31 ±0.03 m w.e./yr; Table 3). During the climatically favourable period in the 1970s and 1980s Zmuttgletscher yielded a smaller mass loss, even below the Swiss-10 wide average (Figure 11). However, although certain areas on the tongue showed an elevation increase (Figure 9c, d, e), the glacier-wide mass balance stayed negative (-0.27 ±0.18 m w.e./yr between 1977 and 1988). After 1988, the mass balance turned more negative again, even if the lowest areas on the tongue still showed stable or even increasing surface elevation ( Figure 4f). The negative trend intensified and since 2001 Zmuttgletscher's mass balance has become more negative than the Swiss-wide average. 15

25
Velocities in the lower part of the glacier tongue (red in Figure 12) have decreased from 10-20 m/yr in the 1980s to almost zero (<3 m/yr since 2005). The central section (orange in Figure 12) showed values of 30-40 m/yr until the 1980s; afterwards velocities strongly decreased to ~5-10 m/yr today. The highest part of the tongue (blue) has decreased over the same period from 50-60 m/yr to ~15 m/yr. In the periods 1977-1983 and 1983-1988 a velocity increase of close to 50% was observed in all areas, with a slightly delayed signal further down on the lower tongue. It was followed first by a strong velocity decrease 5 in the 1990s that got weaker thereafter. Since the 2000s the lowest part of the tongue (red) as well as the Lower SBG show values around 5 m/yr and smaller, which can be considered as almost stagnant. annually. In contrast to the GLAMOS field measurement time series, our results have been acquired with one consistent method, by the same observer and in the same study. All our data points originate from direct measurements on orthophotos, except 1879, which was extracted from the Siegfried map and contains the interpretation of the field analyst and the cartographer, introducing a -presumably small -potential error. 15 The DTMs used for calculating mass balance and surface elevation changes have uncertainties of ~1-2 m (i.e. ±0.2-0.4 m/yr yr -1 ; Table 1), which are derived as the standard deviation over stable terrain and therefore rather represent the upper bound of the uncertainty (Magnússon et al., 2016). Except for the periods 1879-1946 and 1977-1988, these uncertainties are significantly lower than the glacier changes and thus have thus little impact on the final decadal elevation change rates. High elevation change values outside of the glacier surface are found in very steep terrain and in the immediate forefield, which 20 and in the latter case stem from artefacts on water surfaces (e.g. Figure 10 panel e), or the activities of a gravel extraction mining company (e.g. Figure 10g, h, i), and also do not influence affect the results on glacier elevation change.
The time series of glacier surface flow velocities along the lowest ~4 km is over 50 years long  and includes periods of glacier acceleration and deceleration. Even though it covers less than half of the total investigation period, it is one of the longest time series of a debris-covered glacier, similar comparable to the velocity series from to the series by Capt et 25 al. (2016), who extracted velocities for debris-covered Glacier de Tsarmine from 1967 to 2005 (Capt et al., 2016). Even longer time series using single boulder measurements were performed and collected by Kirkbride (1995) and Haritashya et al. (2015) for insights into the evolution of Tasman Glacier since the 19 th century. Despite difficulties in tracking boulders

27
Although it was difficult to retrieve results for boulder tracking in more active areas on Zmuttgletscher due to large displacements (long periods between image dates) and crevasses (e.g. massive increase of in crevasses between 1977 and 1988in 1980s), several representative data points could be found extracted per glacier sector sub-section and time period, and with this information we couldwhich allowed to draw a clear picture of the dynamic evolution of the entire glacier tongue since 1961. Unfortunately, the time period 1946-1961 yielded too few data points (due to image quality and snow cover) and 5 was thushad to be excluded from the time series. In the lowest areas the velocity field derived from the orthophotos 2016-2017 orthophotos (Fig. S14) shows slightly higher values than the interferometer (e.g. 4 ±0.25 m vs. 1 ±0.57 m/yr yr -1 ) and a reduced extent in low flow, , indicating a smaller quasi-stagnant glacier partbut confirms the stagnation of the lower tongue.
This difference can be attributed to the short duration and time of observation in late summer, when flow velocities are likely below the annual average due to the establishment larger uncertainties of the radar (which only records displacements in 10 view direction) as well as to the season of measurements: especially in lowest areas flow velocities might be below annual average in August compared to late spring or early summer, due to the formation of an effective basal englacial and subglacial drainage system (e.g. Nienow et al., 1998).
Debris cover was manually mapped based on topographic maps and orthophotos. A continuous debris cover has reflection characteristics different from ice and can be distinguished visually. The transition from debris-free to debris-covered is 15 continuous, and regions of sparse debris cover -which might even lead to increased ablation -were not included. We mapped only coherent areas of continuous debris cover, thus avoiding to include mixed pixels that indicate a sparse debris cover. Because of varying properties of the orthophotos (RGB vs. panchromatic, resolution, contrast, texture) an automated method would not deliver equally satisfying results and was not applied. The uncertainty of the mapped debris extent is equal to that of the glacier and ice-cliff mapping at, i.e. ±½ pixel. The debris cover of For 1961, the debris cover extent 20 might be slightly underestimated as is less certain, because the orthophoto does not fully cover the ablation area and a is additionally thin snow-covered was present in higher areas. This could partly explain the slight reduction in debris extent from Thus, the area might be somewhat underestimated, which would partly explain the decreasing debris-covered area between 1946 (3.981 ±0.01 km²) andto 1961 (3.6779 ±0.01 km²)., nevertheless we think that the observed trend is correct. Mapping of exposed-ice areas is commonly done manually or by automatic detection including manual correction (e.g. 25 Kraaijenbrink et al., 2016Ragettli et al., 2016Brun et al., 2016). The automaticThe detection of ice cliffs using only topographic and geometric variables as applied by others (e.g. Brun et al., 2016;Kraaijenbrink et al., 2016) showed in our case no satisfying results. , Many steep areas (>40°) were incorrectly identified as debris-covered, and slope values extracted over ice cliffs were often below those typically observed in-situ elsewhere (e.g. Reid and Brock, 2014). e.g. because many steep areas (>40°) are still debris-covered and many cliff areas show up less steep in the DTMs. Our object-based mapping 30 approach is efficient to map larger numbers of cliffs for several datesortho-photos;. However, the polygons are based on a somewhat subjective choice of separation parameters and the manual correction of the polygons is time-consuming. however, the largest time share was still spent on manual selection of automatically delineated polygons, and there is a subjective choice of separation parameters and the selection of polygons. Depending on these aspects and the When considering image quality (texture, contrast, snow cover, cast shadow) the results include both omission and commission 35 errors, which also lie in the range of ± 1 / 2 pixel along the polygon margins.

6.2
Explanation of glacier evolution

Chronological evolution and interactions
The Zmuttgletscher has shown a similar behaviour, pointing out that its overall magnitude of mass changes has 5 foremost been governed by climatic changes, rather than the debris cover.

Chronological evolution and process interactions
On the basis of their dynamic behaviour we distinguish four distinct phases over the entire observation period since 1879. In the first phase, between 1879 and 1977, as a response to the atmospheric warming in the mid-19 th century, the mass balance of Zmuttgletscher turned negative (Figure 3), first only slightly (1879-1946: -0.09±0.12 m.w.e. yr -1 ), then more strongly 10 (1946-1977: -0.67±0.16 m.w.e. yr -1 ), leading to glacier thinning and retreat. At the same time,

15
As an effect of the warming atmosphere since the mid-19 th century, the mass balance of Zmuttgletscher turned negative, first only slightly (1879-1946: -0.09 ±0.12 m.w.e./yr), then more strongly after 1946 (-0.67 ±0.16 m.w.e./yr; Figure 11), leading to glacier retreat and shrinking. On the other hand, thethe debris-covered area constantly increased from 2.758 ±0.2 km² in in 1859 to 3.678 ±0.01 km² in in 1961. The temperature rise isRising air temperatures are probably the main reason for the strong increasinge of debris extent as they -covered area: higher temperatures lead to a rise of the rain-snow transition 20 altitude as well as to more melthigher ice ablation during summer and autumn. Therefore, debris emergence can happen faster accelerates and moves and further up-glacier (Stokes et al., 2007). Decreasing ice flow further also supports the development of a continuous debris cover of substantial thickness as: the emerging debris is evacuated more slowly, thusand the debris has more time to melt out of the ice and can thicken (Kirkbride and Deline, 2013). In the mid-20 th century, a continuous debris cover existed on the lower part of the Zmuttgletscher tongue, likely already responsible for a flattening of 25 the mass balance gradient. Exposed-ice areasIce cliffs and surface flowsupraglacial meltwater channels already existed in the lowest, debris-covered part of Zmuttgletscher already in 1930 ( Figure S137) and thereafter increased in number and area until 1977. These surface features are partly responsible for increased thinning elevation change rates that are higher in these areas close to the terminus than compared to further up-glacier.
In a second phase between 1977 and 1988, Subsequently, aatmospheric cooling and increased precipitation brought a shift 30 towards less negative mass balance from 1977-1988 (-0.27 ±0.18 m w.e. yr -1 , Figure 3/yr). This is reflected in a partial a slight thickening of the tongue (Figure 5Figure 10d, e) as well as in increased ice flow velocities (Figure 11Figure 8, Figure   12b), leading also to enhanced crevassing in the regions of extensional ice flow. ). At the same time, the upper boundary of debris cover in these areas migrated slightly down-glacier (Figure 12a), but because the debris cover continued to expand in Field Code Changed other areas, the total debris-covered area was roughly stable during this period.The total debris-covered area during this period was stable but its upper margin was even moved down-glacier in the central glacier parts where flow velocities increased ( Figure 15). The upper margin boundary of the zone with most exposed-ice areasice cliffs was was also moved down-glacier during this period and the total area of such ice cliffs dropped by a factor of 2 (Table 3Table 3). The reaction adjustment of the of glaciers length and area was less direct, but the retreat and area loss gradually slowed downslowly 5 attenuated.
In the third phase between 1988 and 2001 the After 1988 mass balances turned became more negative again -along with an abrupt rise in air temperature (Figure 3a, b) -more negative, and flow velocity reductions were pervasive across the glacier ies quickly dropped all over the tongue, first higher up, then also in the lowest regions (Figure 8, Figure 12 and b). At Thethe same time, the debris-covered area again strongly substantially increased -especially during the 1990s -while the glacier 10 area was stable and the central part of the terminus even advanced for some metres, leaving behind a small moraine after starting to retreat in 2001 (Figure S 158). The delayed response to effect of the preceding period of more positivewith less negative mass balance is, however, still detectable in the slight elevation gain right at the terminus visible in the elevation change patterns of the tongue until in the period 1995-2001, when the rate of change was still close to zero at the terminus, whereas it was already highly negative on the rest of the tongue (Figure 10f) . 15 In the most recent period (2001 to 2017), the glacier has continued to develop an ever more negative mass balance During the last 1. 5 decades (2001-2017) mass balances have continued on a strongly negative level (approx. ~0.8-0.981 m w.e. yr -1 /yr). Thinning has been most pronounced in debris-free areas -especially on the tongue of STG -and in the ice-cliff zone, andand thinning rates have been constantly remained high (> 1-2 m /yyr -1 r) all overacross the tongue. Towards the terminus these thinning rates are, however, not comparable to values typically observed on debris-free glacier tongues at similar 20 elevations (e.g Findelgletscher: >5 m yr -1 , Sold et al., 2016), which might explain the relatively small Thinning has been most pronounced in debris-free areas -especially on the tongue of STG -and in areas of exposed ice within the debris-

Unprecedented debris cover increase
Zmuttgletscher's debris extent has increased by ~19% points since the end of the LIA without evidence of particularly large rock falls. Here the headwalls surrounding the accumulation area are the most relevant input sources, but lateral moraines 10 also have a minor influence in marginal areas (van Woerkom et al., 2018). Increases in debris-cover extent (for periods from several years to decades) have been shown for other Alpine glaciers (Deline, 2005;Kellerer-Pirklbauer et al., 2008), the Caucasus (Stokes et al., 2007), and the Himalaya (Bolch et al., 2008a;Bhambri et al., 2011), but have so far not been quantified for a century-long period. These studies found debris-cover extent increases of 2-10% points for glaciers with absolute debris-cover extents of 2-42%, which is well below the rates we found for Zmuttgletscher, even on a decadal scale. 15 At Zmuttgletscher we could also identify strong temporal variations in the rate of debris-cover extent change (e.g. stable extent in 1970s and 1980s and high rates of change of >4% per decade in the 1990s and early 2000s), which coincide with distinct changes in climate and related dynamic responses. This strong link to climate may also help to explain the strong debris-cover extent increase on Zmuttgletscher, which we suggest is the result of (i) a combined effect of elevated temperatures and decreasing velocities, (ii) the large debris-free ablation area that still existed a few decades ago, and likely 20 also the (iii) the comparatively small glacier size, i.e. the shorter response time.

5.4
The influence of debris cover on geometry To understand the potential long-term effects of a debris cover on glacier evolution, regional comparisons provide a useful context. Nevertheless, one has to be careful in comparing length and area changes of glaciers for a given climate history, because of differences in geometry and hence response times (Jóhannesson et al., 1989). The retreat of Zmuttgletscher is 25 relatively modest compared to that of other large Swiss glaciers ( Figure 4) and it has shown little terminus fluctuation, even during the climatically favourable period in the 1970s and 1980s. Other glaciers with similarly subdued fluctuations are Unteraargletscher and Glacier de Zinal, which are also debris-covered in their lower reaches. Unlike smaller and debris-free accompanied by a sharp terminus advance in the 1980s (+41 m yr -1 ), followed by a strong retreat since 1985 (-44 m yr -1 ).
During these recent periods of strongly negative mass balance, the retreat rate at Zmuttgletscher (-7.2±0.01 m yr -1 ) was lower than that of all other glaciers in Figure 4.
Area changes from the beginning of the 1970s until ~2012 reveal that Zmuttgletscher has lost 7% of its area (-0.17% yr -1 ), which is comparatively high (Table S15). In comparison, Findelgletscher has lost more than 11% of its area (-0.32% yr -1 ) 10 during a similar period. Thus, in terms of length and area changes, Zmuttgletscher shows a distinctly subdued and delayed response to climatic perturbations compared to debris-free glaciers in the region.
The relatively uniform and reduced thinning over much of Zmuttgletscher's tongue can be attributed to the debris cover and is observed at several debris-covered glaciers in Switzerland. This suppression of thinning by debris cover on the lower ablation area is best illustrated when plotting the elevation change for one period for 11 larger Swiss glaciers along ten 15 vertical elevation bins of equal size ( Figure 13). For all debris-covered glaciers, including Zmuttgletscher, the thinning becomes almost independent of elevation for the lower elevation bins, whereas for debris-free glaciers, it continues to increase towards the terminus. The non-linear relationship between thinning and elevation causes the slower adaptation of debris-covered glacier length and area and hence more extended tongues. A similar flattening or even reversal of the elevation change gradient towards debris-covered glacier termini has also been observed in other regions, e.g. in the central 20 Himalaya (Inoue, 1977;Bolch et al., 2008a;Benn et al., 2012;Ye et al., 2015;Ragettli et al., 2016) and the Tien Shan (Pieczonka and Bolch, 2015). Aletschgletscher (yellow line) is actually also slightly debris-covered.

Climate-driven glacier mass balance and flow dynamics
The evolution of glacier-wide mass balance of Zmuttgletscher since 1879 are comparable to centennial trends on other debris-free glaciers in the Swiss Alps, which are closely related to variations in climate (Figure 3; Schmidli, 2000;Bauder et al., 2007;Huss et al., 2010;Hirschi et al., 2013;World Glacier Monitoring Service, 2017). For Zmuttgletscher, the strong 30 connection of mass balance to climate is exemplified by the almost balanced conditions during colder phases before 1920

Formatted: Caption
Formatted: Left and in the 1970/80s as well as by the strongly negative mass balance in the 1940/50s and in the recent decades (Figure 3; Bauder et al., 2007;Escher-Vetter et al., 2009). This clearly demonstrates that the glacier-wide mass balance of Zmuttgletscher has foremost been governed by climatic changes and has not been strongly influenced by debris cover and changes thereof. This is in contrast to the case of Glacier de Miage in Italy, for which Thomson et al. (2000) related the positive mass balance on the tongue over the 20 th century partly to debris cover. The debris cover of Glacier de Miage is, 5 however, thicker than at Zmuttgletscher.
The relatively thin debris thickness at Zmuttgletscher may on the one hand limit the decoupling from climatic influences, and on the other hand still be sufficient to reduce glacier thinning and terminus retreat. This results in an extended and stagnating glacier tongue that increases the area of ablation (relative to a debris-free glacier tongue) which in turn enhances the sensitivity of the glacier-wide mass balance to climatic changes. 10 The observed temporal variations in ice flow and thickness on the main glacier tongue of Zmuttgletscher also closely reflect variations in climate as exemplified by the acceleration and thickening in the period of more positive mass balance in the 1970s and 1980s. Positive mass balances impacting on flow velocities were also observed at debris-covered Glacier de Tsarmine (Capt et al., 2016) and a number of other debris-free glaciers in Switzerland (Glacier du Giétro, Glacier de Corbassière, Mattmark; Bauder, 2017) and Austria (Hintereisferner, Kesselwandferner;Stocker-Waldhuber et al., 2018). 15 Also at Glacier de Miage, a kinematic wave migrating down-glacier was observed between the 1960s and 1980s that led to a small terminus advance in the late 1980s (Thomson et al., 2000). Also the long-term (1960s until today) reduction in flow velocities on the tongue of Zmuttgletscher is in line with observations from all of the above mentioned glaciers.
All these examples show a direct reaction of flow velocites and hence ice flux to climatic changes. This link is consistent with the principle of mass conservation, regardless of debris coverage, but breaks down when the glacier tongue starts to 20 stagnate as a result of the debris cover slowing down the retreat. Such a dynamical stagnation and decoupling has in the last decade also been observed at the tongue of Zmuttgletscher and is characteristic of strongly debris-covered glaciers (e.g. Bolch et al., 2012).
On large debris-covered glaciers, lower flow velocities are often the result of a decreased driving stress due to flat tongues that result from a sustained reversal of the mass balance gradient (Bolch et al., 2008a;Anderson and Anderson, 2016). At 25 Zmuttgletscher, the elevation change gradient has indeed decreased to almost zero during recent decades. Nevertheless, the glacier surface slope has slightly increased over time, therefore we attribute the decrease of dynamic activity to the reduction in ice thickness and, hence, mass flux and climate. This conclusion is coherent with findings by Dehecq et al. (2019) from glaciers in South and Central Asia. 6.3.15.6 The heterogeneity is even increased by the presence of ice cliffs. overall Debris cover has decreased the 30 elevation change gradient over time to almost zero. The role of Ddebris-related surface features on Zmuttgletscher Whilst the presence of a debris mantle may not substantially influence the overall mass loss rate of Zmuttgletscher, we concludefind that the main reason for the observed heterogeneous thinning pattern is the increasingly extensive and also likely thicker debris cover that is heterogeneously distributed over the tongue. The heterogeneity is evenfurther increased by 35 the presence of ice cliffs. Judging from field Field visits and patterns of surface lowering suggest a close association between ice cliff formation and the presence of DEM differencing, the most relevant surface features leading to exposed-ice areas such as ice cliffs -and hence enhanced ablation -on Zmuttgletscher are supraglacial meltwater channels. Surface melt water is often running off superficially over long substantial distances and several water flow channels are abundant exist all over the glacier tongue, outside of crevassed areas. Inside and alongside these channels bare ice often becomes exposed when 40 water washes away the debris or laterally cuts into the ice (Figure 14Error! Reference source not found.a) and the debris slides off the oversteepened, or because locally differential melting steepens the channel walls until the debris slides off. The location of supraglacial meltwater channels on Zmuttgletscher seems closely related to areas of compressional flow -in flat Formatted: Don't keep with next Formatted: Heading 2 and stagnating areas and troughs between ridges -and has been most pronounced in the lowest ~1.5 km, where also most exposed-ice areas (mostly ice cliffs) occurform. This is consistent with the cessation in e up-glacier migration of the ice-cliff boundary upper margin of this zone stopped in 2005, just below a step in the surface small section of steeper surface slope with extensional flow (Figure 12a at profile distance ~3500 m). Exposed-ice areasIce cliffs also develop when en-and subglacial cavities voids -either slowly or rapidly -collapse and the shear planes of the fractured ice get exposed as also 5 observed at Zmuttgletscher (Figure 16(Figure 14b). We suppose suggest that abundant englacial ablation is effectively creating and enlarging such cavities voids (Benn et al., 2012;Immerzeel et al., 2014;; Figure 14Error! Reference source not found.c) on Zmuttgletscher. Supraglacial ponds are known to have a high potential to increase ablation (Sakai et al., 2000;Röhl, 2008;Miles et al., 2016) and are also potential origins of ice cliffs (Figure 14Error! Reference source not found.d), but these features are not prevalent on their relevance on Zmuttgletscher is low due to their 10 small number and area.  Large terminal ice cliffs also exist at the glacier mouth(s)terminus and are responsible for exacerbated retreat a fast retreat in these locations compared to where the terminus is gently sloping and debris-mantledflattens out and the debris stays on the 20 ice. The situation of a terminal cliff at the glacier mouth combined with terminus retreat is also found e.g. at Gangotri glacier (Bhambri et al., 2012;Bhattacharya et al., 2016), whereas some glaciers with a stable terminus (e.g. Khumbu, Miage) do not show such terminal cliffs (Bolch et al., 2008a;Diolaiuti et al., 2009). Remarkably, depressions and features of irregular surface topography near the terminus are were already indicated in the Siegfried map from 1879. Certain non-terminal ice cliffs on Zmuttgletscher have reached >25 m in height and have persisted for over two decades. Consistent with the general 25 literature, tThe consequences of areas of exposed iice cliffs at Zmuttgletscher are: (i) -a more chaotic pattern of surface elevation changes due strong to locally strong ablation, (ii) -a strongdebris redistribution of debris through cliff backwasting and sediment evacuation and deposition by fluvial transportwater (this redistribution process is difficult to determine but can potentially be important because it leads to a more heterogeneous debris cover), (iii) 30 -and stronger surface elevation changes at the lower part of the tongue (especially downstream of a small topographic steps, see Figure 15a at profile distance ~3500 m). As a result, the elevation change is still stronger towards the terminus than , which would not be the case without these exposed-ice areasice cliffs.
The strong expansion of debris cover might suggest a high and increasing importance of ice cliffs. However, we found that Nevertheless, the exposed-ice areasice cliffs are only responsible for approximately 5% to a small extent responsible for the 35 high volume loss on a of glacier-wide scale (in the range of ~5%) volume loss, because (i) their total exposed-ice area is small (Table 3; also compared to other glaciers: ). On the one hand, exposed-ice areas are smaller compared to other glaciers, e.g. in the Himalaya, where much higher contributions of cliff melt to the total ablation have been observed (e.g. Sakai et al., 2000;Juen et al., 2014;Ragettli et al., 2016;Juen et al., 2014).) and (ii) On the other hand, the debris is relatively thinner debris on on Zmuttgletscher. as well as the conservative calculation method (including buffer areas around exposed-ice areas) might partly explain the smaller elevation change contrast between areas with and areas without exposed ice.
On large debris-covered glaciers, smaller flow velocities are often the result of a decreased driving stress due to flat tongues that result from sustained reversal of the mass balance gradient (Bolch et al., 2008a;Anderson and Anderson, 2016). At Zmuttgletscher, the slope has slightly increased over time, therefore we ascribe the decrease of dynamic activity to the 5 reduction in ice thickness and, hence, mass flux and climate, analogue to findings by Dehecq et al. (2019). The main reason for the observed thinning patterns is probably the increasingly extensive and likely also thicker debris cover that is, however, heterogeneously distributed over the tongue (Figure 7). Due to the stagnant terminus, debris can leave the glacier surface only by terminal retreat or fluvial transport. Decreasing flow dynamics and continuously high ablation will lead to further debris melt-out at the surface, which will increase the debris thickness and in turn decrease ablation. However, in the lowest 10 glacier sections this process is disturbed by surface meltwater channels and ice cliffs and related debris redistribution.

5.7
Dynamic interaction between flow velocities, debris cover, and ice cliffs 20 Even though high flow velocities may be expected to subdue the emergence of exposed-ice areasice cliffs, we could not find a clear velocity threshold linked to their occurrence, e.g. to designate areas of potential ice-cliff formation. Interestingly, dHowever, the upper boundary of the ice-cliff zone was moved downstream by increased flow velocities in the 1970s and 1980s (Figure 12).uring the period of more positive surface mass balance in the 1970s and 1980s, the glacier became more dynamic, and the upper boundary of exposed-ice areas moved downstream ( Figure 15). Further, during this period the 25 absolute and relativetotal ice-cliff area of exposed ice was clearly reduced compared to both periods with lower flow speeds before and after (Table 3Table 3)). This suggests a clear link between dynamic state and the occurrence of ice cliffs and would imply expanding exposed-iceice-cliff areas on stagnating tongues, consistent with the general interpretation of observations (e.g. Pellicciotti et al., 2015).

Formatted: Heading 2
Increased flow velocities also led to a local down-glacier movement of the upper boundary of the debris extent ( Figure 12).
Similarly, Deline (2005) documented a decrease in debris-covered area for Glacier de Miage between ~1890 and ~1920, caused by an increase in surface elevation and flow velocity (Thomson et al., 2000), which confirms that changes in glacier dynamics can directly impact supraglacial debris extent.

6.4
Placing of results in a wider context 5 Length and volume and changes of glaciers for a given climate history are strongly dependent on glacier geometry and, thus, individual response time scales (Jóhannesson et al., 1989). Therefore, one has to be careful in comparing the length change curve of Zmuttgletscher with other glaciers. However, such regional comparisons still provide some useful context. The retreat of Zmuttgletscher is relatively modest compared to that of other large Swiss glaciers ( Figure 3); it has also shown little terminus fluctuation, even during the climatically favourable period in the 1970s and 1980s. Other glaciers with 10 similarly subdued fluctuation are Unteraargletscher and Glacier de Zinal, both of which are also debris-covered in their lower reaches. Nevertheless, their retreat since the 1880s has with 1.7 and 2.5 km (i.e. -21 m/yr and -14 m/yr) been considerably stronger, albeit a similar total length (~10 and ~7 km today, respectively) as Zmuttgletscher (~7 km; GLAMOS, 2018). Unlike smaller and debris-free glaciers, neither Unteraargletscher nor Glacier de Zinal advanced in the 1980s and 1990s (Figure 3). Aletschgletscher is too large and thick to react to short periods of positive mass balance and shows a clear, 15 but increasing long-term retreat trend since the end of the LIA (i.e. -24 m/yr). On the other hand, Findelgletscher, a nearby debris-free glacier of similar size and thickness as Zmuttgletscher and stretching over a similar elevation range, showed a sharp advance in the 1980s (1979)(1980)(1981)(1982)(1983)(1984)(1985): +246 m, i.e. 41 m/yr), followed by a strong retreat since then (1985-2016: -1355in comparison Zmuttgletscher 1983in comparison Zmuttgletscher -2016.2 ±0.01 m/yr). Also during periods of strongly negative mass balance at Zmuttgletscher, its retreat rate has been smaller than that of many other glaciers. Area changes from the beginning 20 of the 1970s until ~2012 for the same sample reveal that Zmuttgletscher has preserved 93% of its area since then (i.e. -0.17 %/yr), which is comparatively high (Figure 3b). During a similar period, Findelgletscher has preserved only ~88.6% of its area (i.e. -0.32%/yr). Thus, in terms of length and area changes, Zmuttgletscher shows a distinctly delayed reaction to climatic perturbations compared to debris-free glaciers in the region.
Comparing Zmuttgletscher to observations from Himalayan glaciers leads to no coherent message about the influence of 25 debris cover on glacier evolution, because of the heterogeneous pattern of glacier length and area evolution. Some Himalayan glaciers have changed more strongly than Zmuttgletscher, e.g. Gangotri glacier (Bhambri et al., 2012;Bhattacharya et al., 2016) and glaciers in North-Western India (Kulkarni et al., 2011; including debris-free glaciers). Most studies, however, yield smaller numbers, which range in the dimension of Zmuttgletscher (e.g. area changes in the Sagarmatha region; Salerno et al., 2008;Bolch et al., 2008a) or even much below, like in Langtang (Kappenberger et al., 30 1993) and Bhutan (Karma et al., 2003). Compared to the evolution of Swiss glaciers, these numbers do not allow a clear distinction between the reaction of debris-covered and debris-free glaciers. To better understand this reaction and ease a direct comparison, longer histories of changes of geometry, climate, and debris cover of individual glaciers need to be considered to account for long response times.
The effect of relatively uniform thinning over much of Zmuttgletscher's tongue can be attributed to the debris cover and is 35 observed at several debris-covered glaciers in Switzerland. By dividing a sample of Swiss glaciers into ten vertically equal elevation bins, the thickness evolution of these glaciers of various sizes, elevations, and elevation ranges can be compared.
For the period 1980-2010 debris-covered glaciers including Zmuttgletscher (Figure 17) showed smaller or constant elevation changes in the lowest elevation bins than higher up. In contrast, debris-free glaciers mostly yielded an increase in elevation change rates until the lowest elevation class. In case of similar geometries of the glacier tongue, such a pattern will lead to a 40 slower adaptation of glacier length and area of the debris-covered glaciers and hence more extended tongues, at least until most of the ice has depleted or the tongue gets separated from the active glacier part at a topographic step with thin ice. A

Comment [n1]: change title
Comment [n2]: other motivation : important to understand whether the obsered evolution is glacier-specific or also visible on other glaciers with somehow similar settings. similar flattening or even reversing of the elevation change gradient over parts of the glacier tongue is a typical feature of debris-covered glaciers, which has been observed also in other regions, e.g. in the central Himalaya (Inoue, 1977;Bolch et al., 2008a;Benn et al., 2012;Ye et al., 2015;Ragettli et al., 2016) and the Tien Shan (Pieczonka and Bolch, 2015). A similar flattening or even reversing of the elevation change gradient over parts of the glacier tongue is a typical feature of debris-covered glaciers, which has been observed also in other regions, e.g. in the central Himalaya (Inoue, 1977;Bolch et 10 al., 2008a;Benn et al., 2012;Ye et al., 2015;Ragettli et al., 2016) and the Tien Shan (Pieczonka and Bolch, 2015).
The mass balance of Zmuttgletscher from 1879-2017 (-0.31 ±0.03 m w.e./yr) is in the range of many other Swiss glacier's mass balance (Bauder et al., 2007). Thus, unlike what Thomson et al. (2000) found for Glacier de Miage in Italy (positive mass balance on the tongue over the 20 th century), the debris cover on Zmuttgletscher does not seem to have a strong influence on the century-long mass change. During the 1970s/1980s, Zmuttgletscher's mass balance was comparable (only 15 slightly negative or even positive) to other glaciers in Alps, and also the strongly negative mass balance of close to -1 m w.e./yr since ~2001 is observed at other large Alpine glaciers (e.g. Findelgletscher, Glacier de Zinal, Glacier d'Argentière; World Glacier Monitoring Service, 2017). Therefore, the century-long and also the decadal reaction of mass change at Zmuttgletscher is comparable to debris-free glaciers in the Alps. In the second half of the 20 th century, the mass balance of Zmuttgletscher (-0.42 ±0.19 m w.e. for 1961-2001) was slightly more negative than most values from the Himalaya, which 20 generally range between -0.1 and -0.35 m w.e. (Dobhal et al., 2008;Bolch et al., 2012;Zhou et al., 2018;Azam et al., 2018).
Also in central Tien Shan glaciers were found to lose mass at similar rates for the period 1976-2009 (ca. -0.27-0.43 m w.e./yr; Pieczonka et al., 2013).. A major reason for the often much lower mass balance values of Zmuttgletscher compared to the Himalaya is presumably the thicker debris layer on Himalayan glaciers exerting a stronger insulation effect (>0.5 m; Patel et al., 2016;Rounce and McKinney, 2014;McCarthy et al., 2017). 25 Stagnant glacier tongues have often been observed at strongly debris-covered glaciers (e.g. Bolch et al., 2012), which is also true for Zmuttgletscher today. The main reason for the slow-down is the decrease in ice flux, rather than the long-term effect of debris cover leading to a flatter glacier tongue. Positive mass balances directly impact flow velocities, which was observed in the 1970s and 1980s at Zmuttgletscher, but also at debris-covered Glacier de Tsarmine (Capt et al., 2016) and a number of other debris-free glaciers in Switzerland (Glacier du Giétro, Glacier de Corbassière, Mattmark; Bauder, 2017) and 30 Austria (Hintereisferner, Kesselwandferner;Stocker-Waldhuber et al., 2018). Also at Glacier de Miage, a kinematic wave migrating down-glacier was observed between the 1960s and 1980s that led to a small terminus advance in the late 1980s (Thomson et al., 2000). Even the long-term (1960s until today) reduction in flow velocities on the tongue of Zmuttgletscher is in line with observations from all of the above mentioned glaciers. All these examples show a direct reaction of flow velocites to climatic changes regardless of debris coverage. However, the extensive debris cover slows down the terminus retreat of Zmuttgletscher and is thus -indirectly -a contributary cause of very low velocities on the lower part of the tongue.
Zmuttgletscher's debris extent has increased by ~19% points since the end of the LIA without evidence of particularly large rock falls. The increase has happened especially by widening and up-glacier expansion of existing longitudinal debris 5 features, along and next to those glacier regions that are partly fed by avalanches. Potentially relevant amounts of debris input originate from lateral moraines that border Zmuttgletscher and its SBG branch in the lower region of the tongue, but their influence is likely constrained to the glacier margin (van Woerkom et al., 2018). Changes -mainly increase -of debriscovered area have been shown for other Alpine glaciers, the Himalaya, and the Caucasus, but have so far not been quantified for a century-long period: Deline (2005) presented a qualitative description of supraglacial debris evolution for several 10 glaciers in the Mont Blanc region since the 18 th century; Kellerer-Pirklbauer et al. (2008) showed that the debris-covered area of Pasterze in Austria increased from 5.4% to 7.3% between 1964 and2000;Stokes et al. (2007) found an increase of 3-6% points at various glaciers in the Russian Caucasus; Gibson et al. (2017) found a slight decrease of debris-covered area on Baltoro glacier from 2001-2012; Bhambri et al. (2011) found an increase of debris-covered area in the upper Gangotri basin but didn't quantify it. Kamp et al. (2011) documented an increase of 10.4 and 7.9 percentage points for two glaciers in 15 Ladakh between 1990 and 2006. The longest time period of quantified debris cover change exists for Khumbu Himal, for which Bolch et al. (2008a) quantified an increase from 39% to 42% between 1962 and 2005, thereby confirming an earlier analysis by Iwata et al. (2000). The values above suggest that the debris cover increase at Zmuttgletscher has been comparatively strong. Main reasons are a combined effect of elevated temperatures and decreasing velocities, the large debris-free ablation area that still existed a few decades ago, and likely also the comparatively small glacier size, i.e. the 20 shorter response time. Similar to Zmuttgletscher in the 1970s/1980s, Deline (2005) documented a decrease in debris-covered area for Glacier de Miage between ~1890 and ~1920, caused by an increase in surface elevation and flow velocity (Thomson et al., 2000). Thus, changes in glacier dynamics can directly impact supraglacial debris extent at the transition between debris-free and debris-covered glacier surface. The debris thicknesses we found (mostly 10-25 cm) were similar to other Alpine glaciers, e.g. Unteraargletscher (Bauder, 2001;Sugiyama, 2003) or large parts of Glacier de Miage (Mihalcea al 25 2008;Foster et al., 2012), but smaller than on many glaciers in the Himalaya and Karakoram with average values of >0.5 m (Nakawo et al., 1986;Patel et al., 2016;Rounce and McKinney, 2014;McCarthy et al., 2017).

Conclusions
This study presents an over >150-year-long reconstruction of glacier the geometry, surface topography, and debris cover from the LIA untilof Zmuttgletscher today,, as well as observations of surface flow velocities over five decades of surface 30 flow velocities. Debris cover extent has more strongly increased than shown elsewhere (from 13.4% of the total glacier surface in 1859 to 31.8% in 2013). We observed the strongest debris-cover extent increase on record, from 12.9±0.9% of the total glacier surface in 1859 to 31.8±0.06% in 2013. Despite the expansive debris cover, we conclude that the evolution of Zmuttgletscher has been mainly driven by climatic changes, both on a decadal and on a century-long scale. The heterogeneous debris thickness distribution reduces ablation from >15% to up to 80% (compared to bare ice) over large parts 35 of the debris-covered area, causing persistent debris-related thinning patterns. Even if though this debris cover is relatively thin (mostly ~15 cm), it has efficiently attenuated and delayed the reaction response to climatic changes in terms of glacier length and area adaptations adjustment compared to other debris-free glaciers in the Alps, both during periods of mass gain and mass loss (overall values of -12.1 ±0.09 m/yr and -0.17 ±0.02%/yr since ~1859, respectively). Hence, elevation change has become (almost) independent from elevation. Despite the expansive debris cover, our observations demonstrate that the 40 evolution of Zmuttgletscher has been mainly driven by climatic changes, both on a decadal and on a centennial time-scale.In contrast, rates of surface elevation change and gGlacier-wide mass balances are throughout the last century comparable to other debris-covered, but also debris-free Swiss glaciers, corroborating a direct link to climate. The prevalence of localised ablative processes such as ice cliff expansion caused surface elevation change in otherwise heavily debris-covered areas, which leads to increased ice loss in the lowest glacier parts, thereby preventing a reversal of the elevation change gradient.
. Nevertheless, debris cover affects the spatial variability of elevation changes: (1) On a tongue-wide scale, the gradient of 5 elevation changes is (almost) constant, as is typical for many debris-covered glaciers. The major reason preventing a more pronounced/inversed gradient is the presence of ice cliffs between the terminus and ~1.5 km up-glacier, leading to enhanced ablation in these areas. Combined with the topographic setting this has led to a slightly increased slope of the glacier tongue over time.
(2) On a more local scale, elevation changes are driven by variations in debris extent and thickness, and the emergence of exposed ice (either originating from channelled surface water flow, or -less important here -from ponding OverallIn general, flow velocities have strongly decreased in the last two decades compared to values before and after a peak in the 1970s and 1980s -similar to other glaciers in the Alps, and the lower 1.5 km of the tongue is almost stagnant today,even though the glacier tongue has slightly steepened over time. Thus, the low flow velocities are mainly due to thinning and 15 reduced ice flux, and only to a minor degree influenced by debris cover. Higher flow velocities between 1977 and 1988 were triggered by a more positive surface mass balance and the related increased mass flux from the accumulation area, which also caused local glacier re-thickening. This increase led to a local slight down-glacier migration of the debris -cover margiboundaryn, followed by a strong debris -cover increase when post-1988 velocities dropped in the 1990s and rising air temperature rise led to increasingly negative mass balance. These observations proof a strong and direct impact of flow 20 dynamics on debris cover extent. Higher flow velocities also moved the upper margin boundary of the ice cliff zone with most ice cliffs down-glacier, temporarily reduced the ice-cliff area, and eventually led to a slight advance in the 1990s.
These findings suggest a clear and direct influence of flow dynamics on debris extent and through the enhancement of icecliff formation in stagnating glacier parts and vice versa. The above described processes and feedbacks are likely valid and relevant for other debris-covered glaciers in the Alps and elsewhere at potentially different rates and magnitudes. In the 25 context of global warming, it is therefore crucial to include these findings in models for glacier projections.

Acknowledgements
This study was funded by the SNSF project #200021_169775 "Understanding and quantifying the transient dynamics and evolution of debris-covered glaciers". We thank M. Fischer and H. Machguth for making the elevation change dataset over all Swiss glaciers available. Thanks to R. Ganarin for conducting and supporting debris thickness measurements. Thanks to 30 the 3G group for the constructive discussions. Thanks to Owen King for proofreading the manuscript.