GERALDINE (Google earth Engine supRaglAciaL Debris INput dEtector) - A new Tool for Identifying and Monitoring Supraglacial Landslide Inputs

. Landslides in glacial environments are high-magnitude, long runout events, believed to be increasing in frequency as a paraglacial response to ice-retreat/thinning, and arguably, due to warming temperatures/degrading permafrost above 10 current glaciers. However, our ability to test these assumptions by quantifying the temporal sequencing of debris inputs over large spatial and temporal extents is limited in areas with glacier ice. Discrete landslide debris inputs, particularly in accumulation areas are rapidly ‘lost’, being reworked by motion and icefalls, and/or covered by snowfall. Although large landslides can be detected and located using their seismic signature, smaller ( M ≤ 5.0) landslides frequently go undetected because their seismic signature is less than the noise floor, particularly supraglacially deposited landslides which feature a 15 “quiet” runout over snow. Here, we present GERALDINE (Google earth Engine supRaglAciaL Debris INput dEtector): a new free-to-use tool leveraging Landsat 4-8 satellite imagery and Google Earth Engine. GERALDINE outputs maps of new supraglacial debris additions within user-defined areas and time ranges, providing a user with a reference map, from which large debris inputs such as supraglacial landslides (> 0.05 km 2 ) can be rapidly identified. We validate the effectiveness of GERALDINE outputs using published supraglacial rock avalanche inventories, then demonstrate its potential by identifying 20 two previously unknown, large (>2 km 2 ), landslide-derived supraglacial debris inputs onto glaciers in the Hayes Range, Alaska, one of which was not detected seismically. GERALDINE is a first step towards a complete global


Introduction 25
There are currently 215,547 glaciers worldwide covering >700,000 km 2 , excluding the Greenland and Antarctic ice sheets (RGI Consortium, 2017). Supraglacial debris covers 4.4% of this glacier area (Scherler et al., 2018) but for many glaciers it plays a critical role in controlling a glaciers response to climate change, due to its influence on surface ablation and mass loss (Benn et al., 2012;Mihalcea et al., 2008aMihalcea et al., , 2008bNicholson and Benn, 2006;Östrem, 1959;Reznichenko et al., 2010).
Extensive debris coverage can alter the hydrological regime of a glacier (Fyffe et al., 2019), with the potential to increase/decrease downstream freshwater availability (Akhtar et al., 2008), and can play a key role in controlling rates of glacier thinning and/or recession, subsequently contributing to sea level rise (Berthier et al., 2010). This supraglacial debris control is thought to be increasingly important with more negative glacier mass balances, with retreating glaciers being increasingly characterised by expanding debris cover extents (Scherler et al., 2011b;Tielidze et al., 2020). The expansion of supraglacial debris cover is due to: (i) glaciological controls such as thrusting and meltout of sub-and en-glacial sediment onto 35 the surface (e.g. Kirkbride & Deline, 2013;Mackay et al., 2014;Wirbel et al., 2018); and, (ii) debris input from surrounding valley walls through bedrock mass movements (Deline et al., 2014;Porter et al., 2010); and, (iii) remobilisation of debris stores, particularly lateral moraines (Van Woerkom et al., 2019). The relative contributions of 'glacially' derived sediment, which may in fact be the re-emergence of glacially modified mass movements (Mackay et al., 2014), as compared to direct subaerial inputs, is highly variable and there is complex coupling between hillslopes and glaciers that varies with relief 40 (Scherler et al., 2011a). However, recent evidence from the Greater Caucasus region (Eurasia) suggests that supraglacially deposited rock avalanches (RAs), attributed to processes associated with climate change, are a key factor in increasing supraglacial debris coverage (Tielidze et al. 2020). Here we focus on the inputs of RAs, high magnitude (> 10 6 m 3 ), perceived low frequency, long runout landslides where there is disparity between current high rates of activity above ice (Allen et al., 2011;Coe et al., 2018) and our ideas of lagged paraglacial slope responses (Ballantyne, 2002;Ballantyne et al., 2014a). 45 Dating of deposits has shown that large RAs are thought to lag ice-free conditions by some thousands of years (Ballantyne et al., 2014b;Pánek et al., 2017). Events cluster in deep glacially eroded troughs and inner gorges relatively low in the landscape (Blöthe et al., 2015). Numerical modelling has shown how considerable rock-mass damage is possible during the first deglaciation cycle (Grämiger et al., 2017); some of the largest inventories highlight a close association with former glacier 50 limits and the source zones of RAs, particularly in the vicinity of glacial breaches (Jarman and Harrison, 2019). However, almost all of our knowledge of past events relies on the presence of in-situ RA deposits. Due to erosional and depositional censuring such deposits are heavily biased to ice-free landscapes where rates of unmodified preservation are higher, although these are still unlikely to constrain true magnitude-frequencies unless rates of geomorphic turn-over are low (Sanhueza-Pino et al., 2011). In glaciated areas supraglacial landslide deposits are rapidly transported away from source areas, removing the 55 simplest diagnostic evidence of a subaerial mass movement process -a linked cavity and debris deposit. Fresh snowfall or wind redistribution can rapidly cover a rock-avalanche deposit that is many kilometres square in area (Dunning et al., 2015).
If this occurs within the accumulation zone the deposit is essentially lost to all surface investigation and non-penetrating remote sensing until eventual re-emergence in the ablation zone, after considerable modification by transport processes. If a RA is emplaced into the ablation zone, censoring may be seasonal, but through time surface transport disrupts initially distinctive 60 emplacement forms (Uhlmann et al., 2013). This supraglacial debris loading represents a glacier input (Jamieson et al., 2015) and can alter glacier mass balance, influence localised melt regimes (Hewitt, 2009;Reznichenko et al., 2011), and glacier velocity (Bhutiyani and Mahto, 2018;Shugar et al., 2012), leading to speed-ups and terminus positions asynchronous with current climatic conditions. Sometimes this leads to moraines that are out of phase with climate, due to the reduction in surface https://doi.org/10.5194/esurf-2020-40 Preprint. Discussion started: 17 June 2020 c Author(s) 2020. CC BY 4.0 License. ablation and surging (or the slowing of a retreat) caused by large landslide inputs (Hewitt, 1999;Reznichenko et al., 2011;65 Shulmeister et al., 2009;Tovar et al., 2008;Vacco et al., 2010).
Currently, the detection of large supraglacial debris inputs -other than through the most common form of ground-based detection, eye-witness reporting -is through the application of optical satellite imagery. This is a labour and previously computationally intensive process, often involving the downloading, pre-processing and manual analysis of large volumes 70 (gigabytes) of satellite imagery. In Alaska, Uhlmann et al. (2013) used ablation zone Landsat mosaics to suggest the frequency of supraglacial debris inputs from large landslides is underestimated, and increasing over time. Seismic monitoring can also be used to detect large debris inputs onto glacier surfaces (Ekström and Stark, 2013) utilising the global seismic network to detect long-period surface waves, characteristic of seismogenic landslides. Seismic methods have identified some of the largest supraglacially deposited RAs in recent times (e.g. Lamplugh glacier RA (Dufresne et al., 2019)) which are compiled in a 75 database (IRIS DMC, 2017), and, when combined with manual analysis of satellite imagery, gives information on duration, momenta, potential energy loss, mass and runout trajectory. However, landslides are challenging to detect using seismic methods and event positional accuracy is limited to a 20 -100 km radius, due to the lack of high frequency waves when compared to earthquakes, further inhibited by the low frequencies and long wavelengths of dominant seismic waves worldwide (Ekström and Stark, 2013). This also results in an inability to detect smaller landslides due to their weak seismic fingerprint 80 (< 5.0 magnitude (M)). Properties of landslides characterised by long runouts onto glaciers are also difficult to extract because their runouts are seismically "quiet", likely due to frictional melting of glacier ice, causing underestimation of event duration and deposit size (Ekström and Stark, 2013). Despite these difficulties, current studies seem to indicate an increase in the rates of rock avalanching onto ice in rapidly deglaciating regions such as Alaska and the Southern Alps of New Zealand, where the majority of recent (aseismic) RAs are associated with glaciers. This increase has been linked to climate warming (Huggel et 85 al., 2012) and potential feedbacks with permafrost degradation (Allen et al., 2009;Coe et al., 2018;Krautblatter et al., 2013).
Here, we present the Google earth Engine supRaglAciaL Debris INput dEtector (GERALDINE): an open-access tool that utilises Google Earth Engine (GEE), and the Landsat data archive encompassing 37 years of optical imagery. The purpose of the tool is to automatically delimit new supraglacial debris additions over wide areas and timescales, which then allows for 90 rapid user-backed verification of inputs from large landslides specifically. GERALDINE is designed to allow quantification of the spatial and temporal underreporting of supraglacial rock avalanches. We describe the methods behind GERALDINE, verify tool outputs against known supraglacial inventories, and, finally demonstrate tool effectiveness by using it to find two new supraglacial rock avalanches, one of which cannot be found in the seismic archives.

Method 95
GERALDINE exploits the capability and large data archive of GEE (Gorelick et al., 2017), with all processing and data held in the cloud, removing the need to download raw data. By default, it utilises tier 1 Landsat imagery (30 m pixel resolution) https://doi.org/10.5194/esurf-2020-40 Preprint. Discussion started: 17 June 2020 c Author(s) 2020. CC BY 4.0 License. that has been converted to top-of-atmosphere (TOA) spectral reflectance (Chander et al., 2009), from 1984 -present, incorporating Landsat 4, 5, 7, and8. GERALDINE also gives the user the following options: (i) to utilise tier 2 imagery; imagery that does not meet the same quality as tier 1, with geodetic accuracy > 12 m root mean square error (Dwyer et al., 100 2018); and (ii) to utilise real time Landsat imagery; imagery that uses preliminary geolocation and where thermal bands require additional processing, before the data is moved to its final imagery tier within 26 days for Landsat 7, and 16 days for Landsat 8. Tier 2 imagery is valuable in regions where tier 1 imagery is limited, e.g. Antarctica where there is a lack of ground control points for imagery geolocation. Real time imagery is useful for rapid identification of landslide locations if a seismic signal has been detected but an exact location has not been identified. Landsat imagery is used in conjunction with the Randolph 105 Glacier Inventory (RGI) version 6.0 (RGI Consortium, 2017). The RGI is a global dataset of glacier outlines excluding those of the Greenland and Antarctic ice sheets, digitised both automatically and manually based on satellite imagery and local topographic maps (Pfeffer et al., 2014). RGI glacier boundaries are delineated from images acquired between 1943 and present day, potentially introducing errors into analysis due to outdated boundaries (Scherler et al., 2018)

(see Supplementary
Information Section 1.0). However, this database represents the best worldwide glacier inventory available and shrinking ice 110 as the dominant global pattern means the tool is occasionally running over ice-free terrain with null results rather than missing potential supraglacial debris inputs. The RGI can be replaced by the user with shapefiles of the Greenland and Antarctic ice sheets, if analysis is required in these regions, or higher resolution (user defined) glacier outlines, if the RGI is deemed insufficient.

Overview of processing flow 115
GERALDINE gathers all Landsat images from the user-specified date range and the year preceding this user-specified date range within the user-specified region of interest (ROI), creating two image collections within GEE (we advise users to define ROIs <5000 km 2 and specify annual date ranges because large ROIs and date ranges can exceed GEE memory capacity). It clips all images to the ROI, applies a cloud mask, then a water mask, before finally delineating supraglacial debris cover from snow and ice. GERALDINE acquires the maximum debris extent from both image collections, creating two maximum debris 120 mosaics, then subtracts these mosaics and clips them to the RGI v6.0 (or user defined area if not using RGI) to output a map.
This map highlights debris within the user-specified time period that was not present in the preceding year, which we term 'new debris additions'. All files can be exported in GeoJSON (Georeferenced JavaScript Object Notation) format for further analysis, including to verify if detections are discrete landslide inputs. An overview of the workflow is presented in Figure 1 and the detail for each step described in Sections 2.1.1-2.1.4. 125

Cloud masking
GERALDINE masks cloud cover using the GEE built-in 'simple cloud score' function (Housman et al. 2018). This pixel-wise cloud probability score allows fast and efficient identification of clouds, suitable for large-scale analysis (Housman et al., 2018) and has been previously applied and well-justified for use in glacial environments (Scherler et al., 2018). A 20% https://doi.org/10.5194/esurf-2020-40 Preprint. Discussion started: 17 June 2020 c Author(s) 2020. CC BY 4.0 License.
threshold is applied to every image, thereby excluding any pixel with a cloud score >20% from the image. We quantitatively 130 evaluated this threshold to ensure optimum tool performance (see Supporting Information Section 2.0). Cloud shadow is not masked as it was found to have a minimal effect on the tool delineating debris from snow/ice whilst greatly increasing processing time.

NDWI mask
Supraglacial streams and/or lakes are present on many glaciers worldwide (for full review see Pitcher and Smith, 2019). 135 GERALDINE includes a Normalised Difference Water Index (NDWI) (McFeeters, 1996) mask to omit these features from debris detection. It utilises the green (0.52-0.6 λ) and near infrared (NIR) (0.76-0.9 λ) bands, the optimum band combination for mountainous regions (Bolch et al., 2011). Other band combinations, such as the modified NDWI that employs the green and shortwave infrared (SWIR) bands (Xu, 2006) as used in Watson et al. (2018), and the blue and NIR band (Huggel et al., 2002), both struggle in the glacial regions we focus on (Chand and Watanabe, 2019;Gardelle et al., 2011). GERALDINE 140 employs a fixed NDWI threshold that was quantitatively evaluated against the glacial lake inventory of Wang et al. (2020), masking values >0.4 as water (see Supplementary Information Section 3.0), similar to previous work on mountain glaciers (Miles et al., 2017). Dynamic thresholding was unsuitable due to processing speed and memory constraints, and, importantly for what the tool is designed to detect, it only offered marginal gains at significant processing and complexity cost.

NDSI 145
The Normalised Difference Snow Index (NDSI) is a ratio calculated using the green (0.52-0.6 λ) and SWIR (1.55-1.75 λ) bands. It helps distinguish snow/ice from other land cover (Hall et al., 1995) and excels at detecting ice where topographic shading is commonplace (Racoviteanu et al., 2008) due to high reflectance in the visible range and strong absorption in the SWIR range. GERALDINE applies the NDSI to all images and a threshold of 0.4 is used to create a binary image of supraglacial debris (<0.4) and snow/ice (>0.4). This threshold has been utilised by studies in the Andes (e.g. Burns and Nolin,150 2014) and Himalaya (e.g. Zhang et al., 2019), but optimum thresholds often vary between 0.5 (Gjermundsen et al., 2011) and 0.2 (Keshri et al., 2009;Kraaijenbrink et al., 2017). We justify our 0.4 value on Scherler et al. (2018) who used this threshold to map and create a global supraglacial debris cover dataset using Landsat 8 images. GERALDINE is in effect standardised with this global supraglacial cover map. We advise users to use this default threshold but if this appears sub-optimum in a user defined region of interest (ROI), the threshold can be fine-tuned in the code (v1.0 line 264 and 274). We utilise NDSI instead 155 of newer band ratio techniques (e.g. Keshri et al., 2009) and more complex algorithms (e.g. Bhardwaj et al., 2015) to ensure transferability between Landsat TM, ETM+ and OLI TIRS sensors as we wish to harness the full temporal archive.

Retrieving maximum debris extent
To attain a maximum debris extent, GERALDINE reduces each image collection to an individual image using a pixel-based approach (Fig. 2). Every binary image (supraglacial debris: 0, snow/ice: 1) in each image collection is stacked, with pixels in 160 the same geographic location stacked sequentially. If any pixel in the temporal image stack is debris, the corresponding pixel in the final mosaic will be a debris pixel, creating a maximum debris extent mosaic. GERALDINE is therefore debris biased due to this processing step. Calculated maximum debris extent mosaics for both the user-defined time period and previous year are differenced, the output being new debris additions. Both the previous year maximum debris extent, and new debris addition mosaics, are displayed for user analysis within the GEE interactive development environment, and easily exportable to Google 165 Drive (included as part of sign-up to Google Earth Engine).

Validation
A bipartite validation was undertaken to assess the effectiveness of GERALDINE outputs for allowing a user to rapidly identify supraglacially deposited RAs: a detection validation (i.e. can the user confirm a rock avalanche has occurred from a GERALDINE output?), and an area validation (i.e. how much of the area of the RA has GERALDINE detected?). Although 170 areal detection is not the main purpose of the tool, greater area detection would ultimately help the user in RA identification.
Validation was performed against the supraglacially deposited rock avalanche (RA) databases of Bessette-Kirton and Coe GERALDINE was run for the year of the event using Landsat tier 1 imagery; the new debris vector output file was exported into a GIS and after an initial qualitative step to see if the user would flag the RA from the GERALDINE output, the area of the deposit it detected was calculated. The tool-detected area was compared against an outline of the RA manually-digitised 180 from a Landsat image using the Google Earth Engine Digitisation Tool (GEEDiT) (Lea, 2018). These two steps allow for an assessment of GERALDINEs ability to highlight new debris inputs, and if this changes over the Landsat era.

Validation
GERALDINE outputs allowed user identification of 92 % of all validation RAs. False negatives all pre-date 1991 (Figure 3), 185 and can be explained by a failure of Landsat satellites from imaging the RA deposit, due to reduced (and insufficient in this case) tier 1 Landsat image availability pre-Landsat 7 within the GEE data catalogue, inhibiting GERALDINE from highlighting the RA as new debris. We note that if just one image featured the RA, GERALDINE would highlight the deposit https://doi.org/10.5194/esurf-2020-40 Preprint. Discussion started: 17 June 2020 c Author(s) 2020. CC BY 4.0 License. as new debris due to its bias towards debris detection (see section 2.1.4). A true 100 % detection rate for RA events on glaciers is, however, unlikely, due to some deposits having high snow/ice content or entraining large amounts of snow/ice during 190 events, which can be common for rock avalanches deposited onto glaciers. This high snow/ice content can mask them as snow/ice during NDSI delineation from debris, inhibiting detection. However, events of this kind also pose significant difficulty for user delineation with original optical imagery. GERALDINE works best when a number of images in the image stack represent maximal debris cover in the preceding year, reducing false positives for the timespan of interest i.e. flagging old debris as new debris, due to a lack of old debris exposure in the previous year. This is particularly applicable to small (<0.5 195 km 2 ) glaciers, where the accuracy of medium resolution satellite imagery is lower (Paul et al., 2013). The debris bias of GERALDINE ensures true negative detection is also extremely high, but this high true negative detection is why user verification of new debris outputs is needed, because they are flagged as new debris but display no supraglacial RA characteristics i.e. lobate and elongated (Deline et al., 2014). To a user familiar with glacial and landslide processes, the differences in GERALDINE outputs between true positives/negatives and false positives/negatives are clear when running the 200 tool to find RA inputs.
GERALDINE RA areal accuracy increases over time from 19 % in the Landsat 4/5 era, to 71 % with the current Landsat 7/8 constellation (Fig. 3), with the latter period characterised by increasingly modern sensors with greater spectral and temporal resolution. Low areal accuracy in the Landsat 4/5 era is once again a product of the GEE data catalogue having limited imagery 205 for certain years in glaciated areas, reducing the ability of GERALDINE to detect the entire area of new debris additions. Areal accuracy increases after the failure of Landsat 4 in December 1993, at which point Landsat 5 is the sole data collector of imagery at a frequency of every 16 days. Despite this single functioning satellite, the tool detects all eight validation events and on average 59 % of the deposit areas between 1993 and the activation of Landsat 7 in 1999. The dual Landsat 5/7 constellation increases tool area accuracy further to 69 %. However, a decrease in mean area accuracy is evident after the 210 failure of the Landsat 7 Scan Line Corrector (SLC) in May 2003 (Markham et al., 2004), decreasing tool areal accuracy by 4 %, due to images missing up to 20-25 % of data per image in the stack (Hossain et al., 2015). We find that a number of Landsat 7 scenes also feature stripes of no data, pre-dating the SLC failure, and can inaccurately cause 'stripes' of new debris in tool outputs. The current Landsat 7/8 constellation has the highest accuracy for detecting the area of RAs at 71 %. The smallest new debris addition we used for validation was 0.062 km 2 , of which GERALDINE detected 71 % of the area, so we have 215 confidence in detection greater than 0.05 km 2 , equating to ~56 Landsat pixels. Even with GERALDINE performing well, additional refinement and/or full automation of RA identification would be an interesting, and priority, area for further investigation.
GERALDINE is frequently affected by the RGI dataset causing over/under-estimation of previous year debris extents and new 220 debris additions. For example, at tidewater glaciers that have undergone retreat since their margins were digitised, the tool often detects clean ice and debris at the tongue. This is solely dependent on the presence of ice mélange in imagery, and NDWI https://doi.org/10.5194/esurf-2020-40 Preprint. Discussion started: 17 June 2020 c Author(s) 2020. CC BY 4.0 License. masking struggling to mask fjord water due to its optimization for masking supraglacial ponds (see Supporting Information Section 1.0 and 3.0). In addition, we found an instance where a supraglacial rock avalanche deposit had been misclassified as a nunatak (60°27'23.7"N, 142°33'35.7"W) and therefore this section of the glacier is erroneously missing from the RGI dataset 225 altogether, preventing tool detection, but this is likely a single case. Topographic shading on debris cover can at times cause pixels to be masked from Landsat scenes due to misclassification as supraglacial water and/or cloud; however, if the tool is run over a sufficiently long period, this will not influence new debris detection.

New Supraglacial Input Detection Example
The Hayes Range, Alaska has a history of large supraglacial debris additions (e.g. Jibson et al., 2006), but no events have been 230 documented in the last decade, in contrast to a recent dense cluster in the Glacier Bay area of Alaska (Coe et al., 2018), which formed part of the validation dataset. To test this, we ran GERALDINE for 2018 to highlight new debris additions on glaciers in the Hayes Range (Fig. 4a). The output map highlighted two large RAs deposited onto glaciers between 1 January 2018 and 31 December 2018. These were manually verified and the potential window of event occurrence identified using satellite imagery within GeeDiT (Lea, 2018). The larger of the two RAs is from a slope collapse on the southern flank of Mt Hayes 235 (4216 m) (63°35'11.7"N, 146°42'50.0"W), with emplacement determined between 10 and 25 February 2018 (Fig. 4b). This rock avalanche was also detected using the seismic method (Ekström and Stark, 2013  covered by snow after the event. The second, smaller RA occurred between 4 and 7 July 2018, on an unnamed glacier to the east of Maclaren Glacier (63°20'21.9"N, 146°26'36.1"W) (Fig. 4c). GERALDINE detected 78 % of this 2.01 km 2 supraglacial debris input, which transformed the glacier from 28 % debris covered to 72 % debris covered, and will have important implications for glacier melt regime, velocity and response to atmospheric drivers. Unlike the larger RA from Mt Hayes, this event was not automatically detected using seismic methods (Goran Ekström, personal communication, 2019), suggesting that 245 its seismic signature was lower than the seismic detection limit (M < 5.0) (Ekström and Stark, 2013). Therefore, there is a high potential to detect all events using GERALDINE, and then provide time-location filters to seismic records to retrospectively quantify force histories and precise timings of events not flagged automatically as a landslide.
We note that new large debris inputs are partially highlighted on the Black Rapids Glacier for 2018 (Fig. 4d), but these 'new' 250 additions were actually deposited in 2002 during the Denali earthquake (Jibson et al., 2006;Shugar et al., 2012;Shugar and Clague, 2011). We assign this discrepancy to minimal cloud-free imagery during summer (a time when deposits are uncovered by snow melt), preventing the tool from highlighting their full summer extent, and causing underestimation of the 2017 debris cover. To a human operator, however, it is clear these debris additions are erroneous because 'new' debris is patchy, with 2017 https://doi.org/10.5194/esurf-2020-40 Preprint. Discussion started: 17 June 2020 c Author(s) 2020. CC BY 4.0 License. debris extent and snow/ice preventing detection of a homogeneous deposit. When GERALDINE is run for multiple years 255 sequentially in the area, the user will have already determined the date of these earlier supraglacial landslides.

Tracking new debris transportation
A secondary use of GERALDINE is tracking existing debris. For example, large supraglacially deposited RAs are transported down-glacier, although often the initial emplacement geometry is characteristically deformed and spread due to differential ablation and ice motion (Reznichenko et al., 2011;Uhlmann et al., 2013). The tool can give an indication of this movement, 260 by highlighting new debris at the down-glacier end of the deposit and differencing the distance of this new debris from the previous year's deposit extent. This gives an indication of deposit behaviour, transient residence time, and glacier velocity, which is often unknown at the temporal resolution of Landsat and complex to calculate in high mountain regions .

265
To demonstrate the evolution of a RA through time, we ran GERALDINE for 2012, 2013, and 2014 for the Lituya Mountain RA in Alaska. This RA occurred on 11 June 2012 and was deposited onto a tributary of the John Hopkins glacier (Geertsema, 2012). The upper portion of the deposit was sequestered into the ice after its deposition in 2012, as is common of debris inputs in glacier accumulation areas (Dunning et al., 2015). However, the deposit toe remained visible on the surface, likely because it was below the snow line. We estimate the down-glacier transport velocity of this RA by tracking and measuring the 270 movement of the deposit toe, to measure the displacement of the deposit leading edge. Using this method, estimates of downglacier transportation of the deposit leading edge between 2012 and 2013 are ~575 ± 30 m, and ~328 ± 30 m between 2013 and 2014 (Fig. 5), the latter in agreement with glacier velocity calculated by Burgess et al. (2013) between 2007 and 2010 (250 -350 m a -1 ), and ITS_LIVE velocity from 2013 (300-400 m a -1 ) (Gardner et al., 2018;Gardner et al., 2019). We suggest that the higher RA deposit velocities between 2012 and 2013 are a result of the immediate response of the glacier to reduced 275 ablation rates and, expansion of debris surface coverage as the RA deposit was redistributed off the protected ice-pedestal formed beneath (Reznichenko et al., 2011).

Conclusion
GERALDINE is the first open-access resource that can rapidly highlight new supraglacial debris additions onto clean ice for a user-specified time and location. Using the output maps it produces, it gives an objective starting point from which a user 280 can identify new debris inputs, eliminating the time-intensive process of manually downloading, processing and inspecting numerous satellite images. We demonstrate its effectiveness by verifying it against 48 known, large, supraglacially deposited rock avalanches that occurred in North America, Europe, Asia, and New Zealand. GERALDINE outputs helped identify 92% of all 48 events, with 100% successful identification post-1991 when image quality and availability increases. We showcase how GERALDINE outperforms current methods that assist with user identification of large debris additions, by identifying two new glacial rock avalanches in 2018, in the Hayes Range of Alaska, one of which could not be detected using current methods. Therefore, the frequency of large supraglacial debris inputs is likely historically underestimated. GERALDINE can become part of the repertoire of tools that enable glacial rock avalanches/landslides to be identified in the past, present, and future. It will improve remote detection and characterisation of these events, to help quantify and evaluate their frequency, spatial distribution and long-term behaviour in a changing climate. 290 Code/data availability GERALDINE code and the validation dataset are available at https://doi.org/10.5281/zenodo.3524414. All other results can be recreated by running GERALDINE in the respective example areas. A guide on how to use GERALDINE is provided in Supplementary Information Section 4.0.

Author responsibilities 295
WS developed the tool and wrote the manuscript. SD made substantial contributions to the conception and functionality of the tool, as well as manuscript editing. NR, SB and JT provided useful guidance on tool functionality and contributed to the manuscript.

Acknowledgments
This work was funded by a Newcastle University Research Excellence Award PhD scholarship to W. Smith. We acknowledge 300 the freely available Landsat datasets provided by the USGS and hosted in the Google Earth Engine data catalogue, and the Randolph Glacier Inventory v6.0 (RGI Consortium, 2017). Hillshade IFSAR DTMs used for figure production were collected as part of the Alaska Mapping Initiative (doi:10.5066/P9C064CO) and are available through the USGS EarthExplorer data portal. We would also like to thank Ryan Dick for his feedback on early versions of GERALDINE.

Competing interests 305
The authors declare that they have no conflict of interests. https://doi.org/10.5194/esurf-2020-40 Preprint. Discussion started: 17 June 2020 c Author(s) 2020. CC BY 4.0 License. Figure 2: Reducer diagram -GEE stacks all images in the collection and undertakes pixel-wise analysis of debris cover, to create a mosaic of maximum debris cover extent. If just one pixel in the image stack is debris, then the corresponding pixel in the maximum debris mosaic will be debris. White pixels represent snow/ice, black pixels represent debris.