Snow avalanche detection and mapping in multitemporal and multiorbital radar images from TerraSAR-X and Sentinel-1

Snow avalanches can endanger people and infrastructure, especially in densely populated mountainous regions. In Switzerland, the public is informed by an avalanche bulletin issued twice a day during winter which is based on weather information and snow and avalanche reports from a network of observers. During bad weather, however, information about avalanches that have occurred can be scarce or even be missing completely. To assess the potential of weather-independent radar satellites, we compared manual and automatic change detection avalanche mapping results from high-resolution TerraSAR-X (TSX) stripmap images and medium-resolution Sentinel-1 (S1) interferometric wideswath images for a study site in central Switzerland. The TSX results were also compared to available mapping results from high-resolution SPOT-6 optical satellite images. We found that avalanche outlines from TSX and S1 agree well with each other. Cutoff thresholds of mapped avalanche areas were found with 500 m2 for TSX and 2000 m2 for S1. S1 provides a much higher spatial and temporal coverage and allows for mapping of the entire Alps at least every 6 d with freely available acquisitions. With costly SPOT-6 images the Alps can even be covered in a single day at meter resolution, at least for clear-sky conditions. For the SPOT-6 and TSX mapping results, we found a fair agreement, but the temporal information from radar change detection allows for a better separation of overlapping avalanches. Still, the total mapped avalanche area differed by at least a factor of 3 because with radar mainly the avalanche deposition zone was detected, whereas the release zone was very visible already in SPOT-6 data. With automatic avalanche mapping we detected around 70 % of manually mapped new avalanches, at least when the number of old avalanches is low. To further improve the radar mapping capabilities, we combined S1 images from multiple orbits and polarizations and obtained a notable enhancement of resolution and speckle reduction such that the obtained mapping results are almost comparable to the single-orbit TSX change detection results. In a multiorbital S1 mosaic covering all of Switzerland, we manually counted 7361 new avalanches which occurred during an extreme avalanche period around 4 January 2018.


Introduction
Snow avalanches frequently threaten people and infrastructure in Switzerland and other mountainous countries. Every winter, dozens of people caught in avalanches suffer serious injuries or even die (Techel et al., 2016) and roads and railways have to be closed during periods of high avalanche danger. To inform about the current avalanche danger levels, ranging from 1 (low) to 5 (very high) on the European Avalanche Hazard Scale (Meister, 1995), the WSL-Institute for Snow and Avalanche Research 25 (SLF) publishes an avalanche bulletin twice a day during winter (SLF, 2018e). The bulletin is written by avalanche experts which analyze weather station data, local snow conditions, detailed weather forecast information and avalanche occurrence reported by a network of in-situ observers. Unfortunately during high avalanche activity low visibility and closed valleys and ski resorts can lead to incomplete or missing avalanche occurrence information. In such situations, as happened in Switzerland in January 2018 and 2019, avalanches can be mapped manually in optical airborne images (Bühler et al., 2009;Eckerstorfer 30 et al., 2016;Korzeniowska et al., 2017) or satellite images which have to be tasked in rapid mapping mode (Scott, 2009;Lato et al., 2012;Bühler et al., 2019). The resulting avalanche outlines can then be used to update avalanche databases which are of great value for hazard mapping and mitigation measure planning (Rudolf-Miklau et al., 2014). As manual mapping is very time-consuming, attempts have been made to automatize avalanche mapping in optical data (Bühler et al., 2009;Lato et al., 2012;Frauenfelder et al., 2015;Korzeniowska et al., 2017). To provide weather-independent observations the project Alpine 35 Avalanche Forecast service (AAF) evaluated terrestrial and spaceborne radar images (Bühler et al., 2014). They concluded that medium to large avalanche events could be mapped using very high resolution radar satellites but with the drawbacks of limited availability and high costs. Nevertheless, for freely available but medium-resolution Sentinel-1 radar images few but promising manual and automatic avalanche mapping studies exist (Vickers et al., 2016;Eckerstorfer et al., 2017;Wesselink et al., 2017;Abermann et al., 2019;Eckerstorfer et al., 2019). 40 To evaluate the applicability of high and medium resolution radar images for avalanche detection in the Swiss Alps we compare 10-meter resolution Sentinel-1 radar images, 3-meter resolution TerraSAR-X radar images, and 1.5 meter resolution SPOT-6 optical images and analyze different methods using multitemporal and multiorbital radar images for two extreme avalanche events which occurred in Switzerland in Jan 2018.
2 Study area and data 45 The study area, shown in Fig. 1, was determined by the spatial and temporal availability of high resolution radar images from the satellite TerraSAR-X (TSX), operated by the German Aerospace Center (DLR). No systematic TSX coverage is available over Switzerland because data are acquired upon request (Werninghaus and Buckreuss, 2010). To cover the two extreme avalanche events around Jan 4 th and 22 nd 2018 ( Fig. 2) with images acquired from identical orbits as good as possible we searched the archive for a sequence of images which limited the study area to the Alps of Uri in central Switzerland. The 50 orbit repeat time defines the dates and limits the revisit time to 11 days for the first event and 22 days for the second one (one acquisition missing). The images were acquired in X-band (9.6 GHz) with the standard TSX stripmap mode (SM) at a nominal single-look complex (slc) resolution of 2.3 × 3.3 m (rg×az). Acquisitions are listed in Table 1. Snow and weather conditions during the two avalanche events are summarized by Bühler et al. (2019). Details are provided by Winkler et al. (2019) andSLF (2018a, b, c, d) (in German). 55 The full TSX scene (black rectangle in the inset in Fig. 1) covers 55×35 km 2 but for the analysis we selected an area of 15.3×8.6 km 2 (red rectangle in Fig. 1) where both the TSX and the validation data  show a very high avalanche activity. The selected area contains steep topography which ranges from 400 -3200 m.a.s.l.. For visualization of results we show in the following only a small subset (blue rectangle in Fig. 1) of the analyzed area.
Radar images of the satellite Sentinel-1 (S1) were analyzed for comparison. S1 images are acquired globally and systemat-60 ically and are free and openly available for download within 24 hours after acquisition (ESA, 2012). Currently, S1 consists of two satellites, S1-A and S1-B, which alternately image central Europe every six days from the same orbit with an slc resolution of 2.7 × 22.5 m (rg×az) in the interferometric wide swath mode (IW). The S1 images, covering 250×170 km 2 , were selected such that they had orbits and acquisition times similar to TSX. Table 1. Satellite data with local acquisition time (CET = UTC+1). Acquisition modes are abbreviated as SM (stripmap), IW (interferometric wide swath), and MS (single pass multi-strip collection). The full list of S1 acquisitions used for the composite of Switzerland is listed in second TSX image was acquired on 2018-01-11, one day before the second S1 image (2018-01-12). On the day in between, the avalanche activity was very low (Fig. 2) and meteorological conditions were relatively stable (SLF, 2018b).
To assess avalanche detection of entire Switzerland, S1 acquisitions were carefully selected from multiple orbits during a 5 day period from before and after the first event (Gray shading in Fig. 2, acquisition details in Table A1).
To analyze the second avalanche event, the SLF ordered optical SPOT-6 images acquired with the single-pass multi-strip 70 collection mode. With this mode the most of the Swiss alps (300×40 km 2 ) could be imaged in a single day (2018-01-24), at a resolution of 1.5 m. These images were visually searched for avalanches by an expert . For comparison we also analyzed TSX data from 2018-02-02, acquired 9 days later.

Radar backscatter physics of avalanches
We detected avalanches based on the radar backscatter signal and their visual appearance (shape). Fig. 3 illustrates a classi-75 fication scheme from (International Commission of Snow and Ice, 1981). The scheme suggests that all avalanche types are composed of three different zones but for some avalanche types (e.g. loose snow avalanches) zones can be difficult to differentiate. The most upslope part is the release area (Fig. 3, blue) with a smooth surface caused by the failure of the weak layer, followed by the zone of transition (purple) with the stauchwall and some deposition caused by the terrain roughness, and finally the tongue-shaped zone of deposition (red) at the bottom which is covered by densely compacted snow granules.

80
Based on snow properties, the different zones show a different radar backscatter signal. In first order scattering physics the total backscatter intensity of a snow pack, σ 0 snow , can be composed of scattering from the snow surface, σ 0 surf , scattering from the snow volume, σ 0 vol , scattering from the ground below the snow pack, σ 0 ground , and scattering from higher order interactions between different structures in the snow pack σ 0 inter . Currently, there exists no specific model tailored to the backscatter properties of snow avalanches (cf. Eckerstorfer and Malnes, 2015, Sect. 5.3), however, general scattering physics from bi-continuous 85 media and rough surfaces can be applied. In that sense, scattering in snow increases with the spatial correlation length of ice Figure 2. The avalanche activity index is the weighted sum of all reported avalanches for Switzerland (Schweizer et al., 2003(Schweizer et al., , 1998. Dry snow avalanches which started high up but were slowed down at medium altitude by wet snow are indicated as "mixed snow" in the legend. Satellite acquisitions dates are indicated by arrows. Images for the multiorbital S1 composite were acquired during the gray shaded periods (see also Table A1). grains (Wiesmann et al., 1998) and also with increased surface-and interface roughness and with decreasing incidence angle θ (Leader, 1971;Fung and Eom, 1982;Kendra et al., 1998).
For plain dry snow of few meters depth scattering at the ground usually dominates the signal because microwaves between 90 1 and 10 GHz are weakly scattered at the snow surface and within the snow volume and penetrate therefore the snow pack to the ground (Xu et al., 2012;Cumming, 1952;Rignot et al., 2001), see also conclusion and simulations in Leinss et al. (2015).
For dry snow the ground roughness determines the backscatter signal but for smooth ground mainly forward scattering (away from the sensor) occurs. For deeper snow or higher frequencies the signal can be dominated by volume scattering (Watte and MacDonald, 1970).

95
In contrast to plain dry snow, snow is deeper and denser in the deposition zone where the surfaces of the avalanche debris can be very rough. Because of the higher dielectric contrast due to the higher permittivity (Matzler, 1996), the contribution of σ 0 vol and σ 0 surf to the total backscatter intensity increases. Both the rough surface and the debris volume scatter radiation more omnidirectional (diffuse scattering) compared to an undisturbed snow pack over smooth ground (specular scattering).
For plain wet snow, however, the incoming radar waves are weakly backscattered at the air-snow interface because most 100 radiation is lost by absorption (Tiuri et al., 1984;Cumming, 1952) and also by forward scattering described by Fresnel coefficients.
As the volume and ground contribution is negligible for wet snow avalanche debris, the dominant backscatter signal results from omnidirectional scattering at the increased surface roughness in the deposition zone of avalanches (cf. Eckerstorfer and Malnes, 2015, Sect. 5.3).

105
Based on the above scattering physics, the zone of origin is very difficult to detect in radar images because the weakly scattering snow volume is reduced without major changes in the surface roughness. The zone of transition should be only sometimes visible, depending on the deposition of avalanche debris. Therefore, mostly the deposition zone can be detected by a brighter backscatter signal and the mostly elongated, tongue shaped geometry.
To obtain a high backscatter contrast with respect to the avalanche surrounding the local incidence angle θ should be far away from zero (i.e. away from layover) to avoid the intense specular backreflection from smooth surfaces. Therefore, the visibility of avalanches in radar images should be much better for slopes facing off the radar. These slopes are also imaged with a higher ground-range resolution δ sr / cos θ which can be close to the slant-range resolution δ sr .

115
All radar products were downloaded in the single look complex (SLC) format. The data were preprocessed with the ESA SNAP Sentinel-1 toolbox and also with the GAMMA software for comparison. The workflow using GAMMA was implemented with Nextflow (Di Tommaso et al., 2017) to speed up execution and code development and to ensure a reproducible analysis. Preprocessing consists of coregistration, multilooking for reduction of radar speckle (TSX: 6×5 px, S1: 4×1 px), orthorectification, and generation of radar shadow and layover masks. The SNAP workflow for S1 images is shown in Fig. A1. We did not apply 120 any radiometric terrain correction as the visible topography helps to identify the avalanche path direction.
For orthorectification we used the Swiss elevation model SwissAlti3D (2013) downsampled from 2 m to 30 m resolution.
We noticed, however, that despite of using the same DEM and output resolution, sharp topographic features seem to be better orthorectified with the GAMMA software which might use a more precise spatial interpolation. The radar images were orthorectified to a resolution of 5x5 m (TSX) and 15x15 m (S1) and the backscatter signal in dB was saved to geotiff files. The 125 exact radiometric normalization is irrelevant, because we did not apply any radiometric terrain correction (Small, 2011) and different ellipsoidal corrections (σ 0 E , γ 0 E ) differ only by almost constant factors. Since the TSX data was acquired with a single polarization (the co-polar channel HH) we also used only the co-polar channel (VV) of the two available polarizations of S1 to obtain a fair comparison. For the multiorbital composites, we used both polarizations of S1 (VV, VH).

130
Although avalanches could be manually detected in single radar images they are difficult to analyze with automatic methods. As radar systems carry their own illumination system the backscatter signal is primarily determined by topography and land cover type. It is therefore common practise to analyze change detection images to separate sudden backscatter changes from stable topographic and land cover features (Wiesmann et al., 2001;Eckerstorfer and Malnes, 2015). To correct for large-scale backscatter changes due to wet snow a 500 m highpass filter was applied to the backscatter difference between two consecutive images.

135
Examples for TSX and S1 are shown in Figs. 4a and 4b. To create the images, the backscatter intensities in dB were normalized by clipping the lower and upper 1%. Consecutive images were then stored in the channels [R, G, B] = [img2, img1, img1] so that backscatter changes are well visible by the red-cyan contrast. From the images (TSX and S1) avalanche outlines were drawn manually.
(d) Multiorbital S1 change detection image as in (c) but with nonlocal mean filter applied. The change detection images allow for a temporal classification of avalanches into three classes (new, old, unsure). New avalanches appear red because of increased backscattering and are therefore assumed to have occurred between the first and the second acquisition. Old avalanches, with a decreasing backscatter signal, appear blue are therefore assumed to have occurred before the first acquisition. Bright features with unchanged backscatter intensity appear almost white and are classified as unsure if they look like avalanches.

145
The free and systematic availability of S1 radar images and the short revisit period of six days allow for creation of an RGB composite change detection image covering entire Switzerland. Therefore, 12 images, acquired between 2017-12-28 and 2018-01-01 from different orbits, were combined into an image before the first avalanche event (Jan 4 th ). Another 12 images, acquired between 2018-01-09 and 2018-01-12 with an identical imaging geometry, were used for the post-avalanche event image. The images (listed in Table A1) were preprocessed according to Sect. 4.1. To reduce radar speckle we averaged both polarizations 150 and weighted the cross-pol channel (VH) by the ratio a of the co-and cross-pol backscatter intensities averaged over the acquisition footprint: Then, the weighted mean was converted to dB and scenes from different ascending and descending orbits were averaged.
Thereby, a relatively homogeneous bright image is generated where layover areas lighten up the relatively dark slopes facing 155 away from the radar without screening too much of the contained details (Fig. 4c). To further reduce noise but to preserve edges in the mosaic images, we applied a non-local mean filter (Jin et al., 2011;Condat, 2010). The filtered image is shown in Fig. 4d.

Relative brightness of snow avalanches
To analyze the brightness of avalanches relative to their surrounding, we calculated the ratio of the mean backscatter signal of 160 an avalanche area and its surrounding area. Therefore, a manually generated avalanche mask was dilated once by 9 and once by 18 pixels. The difference of the two masks defines the surrounding. For the avalanche mask, the visual avalanche mask was eroded by 3 pixels to reduce manual contouring errors. To obtain statistically significant results we calculated the backscatter ratios only for avalanches and surrounding areas larger than 100 pixels.

Automated avalanche detection 165
As manual avalanche mapping is time consuming, a reliable automation of this process would make the mapping data quickly available for further application. Therefore, different attempts have been made to automatically detect avalanches mainly on the two satellite platforms S1 (Vickers et al., 2016;Wesselink et al., 2017;Abermann et al., 2019;Eckerstorfer et al., 2019), and Radarsat-2 (Hamar et al., 2016;Wesselink et al., 2017). The general workflow in these papers is quite similar to ours.
All methods are based on two-image change detection, application of various masks (layover, shadow, water bodies, forest), thresholding and filtering of extracted avalanche properties.
In addition to a shadow and layover mask, we applied a slope dependent mask to limit the detection to potential avalanche deposition zones for which we expect the strongest backscatter change. By definition, friction is larger in the deposition zone than the downhill-slope force. Therefore, slopes with an inclination larger than 35 • , which typically occur in the zone of origin, are masked out (Bühler et al., 2009).

175
For noise reduction but to preserve avalanche edges, a 5×5 px median filter was applied to the backscatter difference images in dB. As avalanches should have a well defined edge, an edge mask was generated by applying a Sobel filter with a 5×5 kernel to the median filtered difference image.
In the median filtered difference image, from all pixels brighter than a threshold of 4 dB, the brightest 5% were considered as the mask of potential avalanches. The threshold was determined empirically based on TSX data but other authors also used 180 thresholds of 4-6 dB Karbou et al., 2018;Vickers et al., 2016). To remove isolated bright pixels from the mask, we determined around each continuous area an ellipse and removed areas with a major axis shorter than 15 pixels (TSX: 45 m, S1:225 m). Additionally, only potential avalanches for which more than 10 pixels intersect with the edge mask were considered for the final avalanche mask.

185
None of the mapping results obtained from TSX, S1, or SPOT-6 can be considered as real ground truth and different avalanches or avalanche shapes were detected with the different methods and satellites. Also, ambiguous relations can exist when a single large avalanche in one mapping result appears as multiple smaller avalanches in another mapping result. This makes the evaluation of binary classifies (e.g. probability of detection or false discovery rate) difficult or even impossible. We refrained from using a pixel-to-pixel comparison which would have demanded a manual mapping precision on the pixel level which 190 contradicts the subjective mapping by an experienced expert who sometimes estimates an avalanche outline from discontinuous avalanche patches.
As a remedy we compare results from two data sets A and B by reciprocal counting of avalanches which overlap in both data sets (considered as "found") and avalanches which do not overlap (considered as "not found"). These numbers differ depending on the direction in which the comparison is done (A→B or B→A). Depending which data sets is considered as ground truth, 195 avalanches which were "not found" can be either regarded as false negative alarms (missed) or as false positive alarms (false alarm).
For conciseness we abbreviate the RGB change detection images by acquisition month and day (mm-dd / mm-dd). In the change detection image tsx(12-31 / 01-11), covering the first avalanche period, a total of 267 avalanches were manually detected in the study area (red polygon in Fig. 1). As detailed in Table 2, 164 avalanches were classified as new and 19 were classified as old avalanches. For 84 avalanches a clear assignment to new or old was not possible. Therefore, we assigned them to the class unsure. For example, in the upper part of Fig. 4a arrows indicate two large new avalanches which completely covered two small old avalanches. Therefore, their backscatter signal did not change and these old avalanches were classified 205 as unsure (though they could be classified as old using spatial context information).
In the change detection image tsx(01-11 / 02-02), covering the second avalanche period, a total of 351 avalanches were detected, composed of 170 new avalanches, 35 old ones and 146 unsure cases. Most of these unsure avalanches were actually classified as new after the first avalanche period but overrun by new avalanches during the second avalanche period (compare Fig. A2 with Fig. A3). Therefore, the number of old avalanches seems to remains low.

TSX compared to optical SPOT-6
The SPOT-6 images were acquired immediately after the second avalanche event in the morning of Jan 24 th . Avalanche were mapped by E. Hafner in SPOT-6 images . They found that only 24% of outlines were clearly visible; 76% of the avalanches outlines were estimated between partially visible release and deposit areas. In the study area, the SPOT-6 avalanches did not contain any age information but the authors conclude that 20 − 45% of avalanches were already released 215 many of them are actually old avalanches.
Due to the 11 day revisit time of TSX, the next available TSX image was acquired 9 days after the second event in the evening of Feb 2 nd . During this 9 days surface melt occurred followed by about 20 cm of new snow on Feb 1 st . Due to changing snow properties the contrast between avalanches and the surrounding snow has very likely decreased. Without knowledge of the SPOT-6 mapping results, avalanche were mapped in the TSX images independently by the second author of this work. The 220 outlines differed significantly, however, most likely because different features (avalanche origin, path, deposit zone) are visible in optical and radar images. Therefore we decided for a feature-based comparison, i.e. overlapping polygon are considered as detected in both data sets. Avalanches split up into discontinuous polygons were counted separately, even if all polygons overlap with one single large polygon in the other data set ( Sect. 4.6).   Vice versa, 44% (125/286) of the optically detected avalanches were also found in the TSX change detection image (Table 3b) but more than half of the optically detected avalanches were not found. 20% (57/286) could not be found because they are located in the radar shadow and 36% (104/286) had a too low backscatter contrast to be visible with radar. We did not found 230 significant differences in area for the lower detection limits: for both, TSX and SPOT-6 the smallest detectable avalanches had an area of 500 m 2 (Sect. 5.6).
Using the temporal information from radar change detection the 125 avalanches detected with SPOT-6 but also with TSX (Table 3b) can be further classified into 27 new, 38 unsure, and 6 old avalanches. The remaining 54 avalanches could not be unambiguously classified, because they cover areas differentiated into multiple different classes by radar whereas such a 235 temporal classification is difficult with single SPOT6 images . Figure 5 shows a subset of the SPOT-6 images and visualizes the manually mapped avalanches. Especially in the lower part of the image, in the cast shadow, many small radar-detected avalanches (red) were not found in the optical analysis (blue).
With radar, avalanches could generally not be detected in the radar shadow or layover (added with black) but also many other avalanche were missed by radar.

TSX compared to S1 change detection
To asses the added value of high resolution TSX images we compared them to medium resolution S1 images. We chose the first avalanche period to simplify counting because of less overlapping old and new avalanches. In the S1 change detection image S1(12-31 / 01-12) a total of 89 new, 13 unsure, and 16 old avalanches were found. The S1 image shows a significantly lower resolution than TSX (Fig. 4a vs. Fig. 4b) such that small avalanches (yellow in Fig. 6) are more likely to be missed.

245
As detailed in Table 4, from the 89 new avalanches, 83 were also found by TSX. They correspond to 76 new and 7 unsure avalanches; 6 avalanches were not found. Vice versa, two thirds (104/164) of the avalanches found in tsx(12-31 / 01-11) corre- spond to the 83 avalanches also found with S1. One third (60/164) was not found, mostly because they were too small to be detected with S1. We found that the smallest avalanches detected by S1 have an area of around 2000 m 2 (Sect. 5.6). By combining S1 acquisitions from multiple ascending and descending orbits, we minimized areas affected by radar layover (areas with radar shadow appear in layover when imaged with the opposite pass direction). The multiorbital change detection composite is shown in Fig. 7 and covers entire Switzerland during the first avalanche period. In the full, non-local mean filtered 15 m-resolution image, which is available online (Leinss et al., 2019), we manually counted 7361 avalanches but did not draw any polygons. We found that avalanches reaching below the wet snow line were much better visible than avalanches from the 255 dry snow zone. The subset shown in Fig. 4c illustrates the mitigation of layover (in the upper and lower right), the speckle reduction and the enhanced resolution compared to the single orbit S1 image in Fig. 4b. Only areas near radar shadow loose contrast and show a reduced avalanche visibility because the added layover image does not contain useful information.
The comparison of the multiorbital S1 mapping results with the high resolution TSX data is detailed in Table 5. In the study area a total of 136 new avalanches were manually detected in the multiorbital image (S1-MO). Of these, 104 avalanches match 260 with avalanches detected in the corresponding single orbit TSX change detection scene (95 of them with new avalanches, 9 with unsure), whereas 32 avalanches were not found with TSX. 17 of the 32 avalanches could not be detected because they are in the shadow/layover areas of TSX. Vice versa, 110 of 164 TSX avalanches were also detected in the multiorbital S1 composite whereas 54 TSX avalanches were not detected.  thirds of the manually identified avalanches in the same image pair. Here we considered the manually determined avalanche mask as a proxy for the true extend of the deposition zone. We are aware that the significance of such a comparison is limited. Lines indicate the thresholds for the backscatter difference (dashed) and the slope-dependent mask (solid).

Automated avalanche detection
Nevertheless, the advantage of this comparison is that the performance of the detection algorithm is directly compared to the results of a human avalanche mapping expert.
For the first image pair tsx(12-31 / 01-11) Table 6a details that 110 of 164 manually mapped new avalanches were also 275 found with the automated detection whereas 54 were not found. As shown in Fig. 9, these "missed" avalanches are often small avalanches which were filtered out by the algorithm. Vice versa, of 138 automatically detected avalanches 21 were not found manually (Table 6b).
When considering the total number (164) of manually detected avalanches in the study area (red polygons in Fig. 1) as truth one can assign avalanches which were also found automatically to true positive (TP = 110), i.e. correctly detected. The remain-280 ing avalanches, which were not automatically detected, are then assigned to false negative (FN = 54), i.e. incorrectly rejected.
With this assumption the probability of detection (POD) and the miss rate or false negative rate (FNR) can be calculated: Further, one can assign automatically detected avalanches which were not manually found to false positives (FP = 21), i.e.
incorrectly detected. When assuming that the number of correctly detected avalanches is given by TP = 110, the false discovery With that one obtains a POD = 67 %, a miss rate FNR = 33 % and a false discovery rate FDR = 16 % for tsx(12-31 / 01-11).
The automated algorithm was also run on the S1 images pair S1(12-31 / 01-12). As detailed in Table 6e, 68 of 89 manually 295 detected new avalanches were also found automatically whereas 21 were not found. Vice versa, of 92 automatically mapped avalanches 72 were also found manually and 20 were not found (Table 6f), resulting in a POD = 76%, a FNR = 24%, and a FDR = 23%.
The higher POD and lower FNR for S1 compared to TSX indicates only, that with S1 the automatic method can detect a larger fraction of the manually detected avalanches. It does not indicate that results obtained from S1 are better compared to 300 TSX data where in total more avalanches were detected. : the cumulative avalanche area Σ i 1 Ai plotted over avalanche size (Ai) reveals that the smallest avalanches size detected by TSX and SPOT-6 is about 500 m 2 , 2000 m 2 for S1, and around 1000 m 2 for the automatic methods. The total cumulative areas differ by an order of magnitude: with radar only bright deposit areas of new avalanches were mapped automatically (1.3 · 10 6 m 2 ), and less-bright areas were added manually (2.5 · 10 6 m 2 ). Summing all classes (new, old, unsure) in TSX images results in 7.5 · 10 6 m 2 which is one third of the cumulative area of SPOT-6 outlines (2.5 · 10 7 m 2 ) which cover also older avalanches and for which outlines of the entire avalanche area (release, path, deposit) were either clearly visible or were at least estimated.

Size distribution of detected avalanches
The size distribution of detected avalanches depends on sensor resolution and also on which features are actually visible by the sensor. For radar sensors it is likely that only the deposit area is mapped, whereas for the SPOT-6 dataset care was taken to map (or at least estimate) the entire avalanche area, including the release area (Sect. 5.2). Because with radar only partial 305 areas were mapped, simple size distributions (Fig. 10a) may appear shifted. To provide more detailed insight we plotted the cumulative area Σ i 1 A i of all avalanches sorted by their apparent area A i (Fig. 10b). The smallest detectable avalanche size can be found in the lower tail of the curves in Fig. 10b: for TSX and SPOT-6 the smallest avalanches have about 500 m 2 , 2000 m 2 for S1, and around 1000 m 2 for the automatic methods.
It may surprise that in the study region the total avalanche area in SPOT-6 images is an order of magnitude larger (2.5·10 7 m 2 , 310 green curve in Fig. 10b) than the total area of manually radar-detected new avalanches from TSX and S1 (red, blue, and orange dots: 2.5·10 6 m 2 ). A factor of three remains when comparing the area of all (new, old, and unsure) avalanches detected by TSX (purple in Fig. 10b) with SPOT-6 which does not contain any age classification. Considering the fact that with radar mainly the deposition zone can be mapped a difference of a factor of three is reasonable.

Radar change detection images
The temporal information from radar change detection makes it possible to differentiate relatively clearly between new and old avalanches, at least for low avalanche activity where old avalanches are rarely overrun by new ones. This can be seen as a major advantage compared to optical images for which temporally dense time series are not reliably available due to weather conditions. The missing temporal information can lead to an overestimation of the avalanche area and  320 report that deposit areas of large avalanches (> 10 000 m 2 ) remain visible for several weeks.
Nevertheless, for strong avalanche activity, the differentiation of overlapping avalanches is difficult even with radar. For example, we found a large number of unsure avalanches for the second analyzed avalanche event (Sect. 5.1) which could be assigned to new avalanches of the first event. For temporal separation, fast repeat times of current radar satellites, like 6 days when combining the two S1 satellites is a major advantage compared to other satellites (TSX: 11 days, Radarsat: 24 days). To 325 differentiate overlapping avalanches a recently developed age-tracking algorithm showed promising results .

Optical mapping vs. radar change detection
Regardless of the advantages of radar change detection, the spatial resolution of optical sensors is better compared to the nominal resolution of radar sensors because the intrinsically coherent SAR imaging method makes radar speckle unavoidable 330 and requires spatial or temporal averaging. Furthermore, the resolution of TSX and S1 is not good enough to recognize flow structures of the avalanche surface which are well visible in the optical SPOT-6 images .
Nevertheless, using TSX change detection we have mapped a similar number of avalanches (316) compared to the results from optical SPOT-6 images (286 avalanches) within the study area. However, the mapped avalanche outlines differ relatively strongly and are sometimes split up into sub-polygons which results in the fact that only 68% of the radar detected avalanches 335 overlap with avalanches found with the optical data or inversely, only 44% of optically detected avalanches were also found by radar. The differences of avalanche outlines could also be partially attributed to the fact that some avalanche outlines were estimated by the (different) persons mapping the avalanches.
The fact that a larger fraction (68% vs. 44%) of radar-detected avalanches matches with optically detected ones results from the better differentiation of adjacent avalanches into multiple classes (new, old, unsure) which were often mapped as one large 340 avalanche with optical data. When multitemporal optical data is available, a temporal differentiation is also possible  which, however, was done for a different region than our analyzed area.
From the analysis of avalanches detected by radar but not found in the optical SPOT-6 image, we found that more than 80% of these avalanches were located in the cast shadow. Similar, in radar images no information is available from the radar shadow and very poor information is available from layover areas, however, only 35% of avalanches not found in the radar images (but 345 in optical) are located in the radar shadow or layover. We think it is an important result that not only radar acquisitions are affected by (radar) shadow but that avalanche mapping using optical data seem also to be deteriorated by the cast shadow from high mountains.
A main difference between SPOT-6 and radar mapping results is that the total avalanche area differed at least by a factor of three (Fig. 10b). We attribute this difference to the fact with SPOT-6 avalanches were mapped more completely (origin, path, 350 deposition zone) than with radar (mainly deposition zone). This has important consequences when comparing avalanches by pixel area rather than by overlap.
Due to unfortunate acquisition timing, the direct comparison of SPOT-6 and TSX data is not ideal: the SPOT-6 images (2018-01-24) was acquired just between the two TSX images (2018-01-11, 2018-02-02) which leaves 9 days where additional avalanches could have occurred, considering about 20 cm of fresh snow on Feb 1 st . Nevertheless, Fig. 2 indicates that the 355 biggest part of avalanches occurred before the SPOT-6 acquisition and only about 5% of avalanches occurred until Feb 2 nd . We confirm this by analyzing a multiorbital S1 change detection image S1(01-24+01-28 / 01-30+02-03) where we did not find any new avalanches in the study area.

TSX compared to S1 change detection
The comparison of TSX and S1 change detection images, both of them acquired for the first avalanche period with almost 360 identical orbits and acquisition times, shows that the S1 satellites are a valuable source of radar images for avalanche mapping.
The size of smallest detectable avalanches for TSX are "medium" avalanches (500 -10 000 m 2 ) with a width of more than 20 m. S1 misses mainly "medium" avalanches smaller than 2 000 m 2 (Fig. 10b). Similar results for S1 with a minimum cutoff of 4 000 m 2 were found by Eckerstorfer et al. (2019).
Still about two thirds of avalanches detected and classified as new with TSX could also be detected with S1 (Sect. 5.3). and S1 mapping results. This is confirmed by Fig. 10b which shows that, despite of a different lower cut off, the total area of radar-mapped new avalanches agrees very well (2.5 · 10 6 m 2 ). Also, the shape of the avalanches masked in S1 data is very similar to the one from TSX (Fig. 6). Therefore, we consider the reduced resolution and separability of avalanches in S1 images to be much less relevant than the superior availability of S1 data.

Multiorbital composite
The combination of radar images acquired with different polarizations and from ascending and descending orbits reduced radar speckle and minimized areas affected by layover. By combining two orbit and (pairwise incoherent) polarizations, areas visible from both orbits were imaged by 4 independent observations. In our case of mapping entire Switzerland for a specific period, this number was even increase to 6 or 8 observations when acquisitions with different incidence angles (from the same 375 orbit direction) overlap. Due to the 4-8 independent observations, spatial multilooking (used for speckle reduction) could be reduced to 4×1 pixels to obtain a radiometric accuracy otherwise only possible with multilooking windows of 8...16×2 pixels.
With this multiorbital averaging method, we estimate that an effective spatial resolution of about 20×20 m was achieved (TSX: about 10×10 m after multilooking). This resolution enhancement can be clearly observed when comparing Fig. 4b with Fig. 4c.
Also, about twice as much medium size avalanches were detected compared to a single S1 image (Fig. 10a). However, because 380 topography was not considered during averaging the resolution can deteriorate in slopes facing off the radar (Sect. 3).
Another drawback of combining acquisitions from multiple dates is that no unique time stamp can be given to the "before"and "after"-aquisition. In the worst case, avalanches loose contrast if they had occurred during the collection period of the set of "before" images. However, in our case, we focused on the extreme avalanche event on 2018-01-04 (Fig. 2) and made sure that the "before"-and "after"-imaging period did not overlap with the main avalanche event. For an operational use combined 385 (asc+desc) acquisitions must be acquired within a time-period as short as possible, i.e. significantly shorter than the orbit revisit time to avoid reduced visibility by averaging out "in-between" avalanches which are only visible in one of the two averaged acquisitions. For S1, ascending and descending acquisitions with only 12 hours time difference should be used, if possible.
Considering a revisit time of 6 days over Europe results in probability of 1:12 to reduce the visibility of averaged avalanches.
In this study we simply averaged the change detection radar images and did not apply any terrain correction. We think that 390 more advanced methods to merge radar images from multiple orbits, for example local resolution weighting (LRW) by Small (2012), should further improve avalanche mapping results. From the comparison to optical data we also found that avalanches can be clearer identified in slopes facing off the radar compared to slopes which are facing towards the radar (but not yet in layover). As detailed in Sect. 3, we think that, because of the more isotropic scattering from the rough avalanche debris surface, steep local radar incidence angles should be used to enhance the local contrast to the surrounding snow. Therefore, 395 slopes facing away from the sensor should be given more weight which is done already implicitly by LRW. Furthermore, in mountainous regions LRW applies already unequal weights for ascending and descending acquisitions which decreases further the probability that avalanche falling inbetween the averaged acquisitions loose their visibility.

Automated avalanche detection
For both, TSX and S1 images the implemented avalanche detection algorithm performs with reasonable results, at least when 400 the number of overlapping avalanches is low. That means that in general a few sparse events are more likely to be detected than overlapping clusters of avalanches.
Compared to the manually detected avalanches (red shading in Fig. 9), the area of automatically detected avalanches (yellow) shows a good agreement. However, the upslope parts of avalanches are often only fractionally detected because of their relatively low brightness. For a weakly visible starting or transition zone a human observer can conclude that it must belong to 405 the below situated avalanche deposit. Also, by choosing a threshold of 4 dB already 18 % of the manually detected avalanches are likely to be missed (Fig. 8a). A dynamic threshold based on backscatter changes in individual image pairs could improve these results . Further, minor parts of manually detected avalanches are located in slopes steeper than 35 • (Fig. 8b) which were masked out by the automatic method.

410
The fact that no real ground truth exists makes a direct comparison of the different methods difficult. However, some methods show a much higher potential to differentiate large connected avalanche patches into multiple smaller ones than other methods.
Therefore we use a reciprocal, two-way comparison of avalanche detection numbers to estimate which of the methods can better differentiate adjacent avalanches.
As a proxy for the enhanced differentiation we define the ratio N A→B /N B→A where N A→B is the number of avalanches from 415 data set A which were also found in the data set B, and inversely, N B→A is the number of avalanches in B which were also found in A. Additionally, we define the ratio N A¬B /N A of avalanches found in A but not found in B relative to all avalanches found in A and analogue N B¬A /N B . The meaning of the last two ratios depends on interpretation and correspond to the false discovery rate (FDR) under the assumption that B is considered as truth or alternatively to the false negative rate (FNR) if A is considered as truth. 420 Table 7 lists the three ratios for different data sets. We interpret these numbers such that a differentiation ratio NA→B NB→A > 1 indicates that set A provides spatially more detailed results than set B. An asymmetry between the last two columns indicates that one method detects more avalanches than the other method.
From the comparison to SPOT-6, derived from Table 3, we infer that TSX change detection allows for a better differentiation of avalanches than single optical images. However, both methods show miss rates (and possibly some false detection) of 32 425 and 56% for avalanches which are not visible by the other method which indicates a certain complementarity of optical and radar images for avalanche detection.
Compared to S1, the higher resolution of TSX allows for a 25% better differentiation and 37% more avalanches were detected (derived from Table 4). Still, the false discovery rate of S1 compared to TSX is quite low (7%).
Interestingly, the avalanche separability of the multiorbital S1 composite, including the NL mean filter, is very comparable 430 to TSX single orbit change detection (1.06) while 33% or 24% of avalanches detected by one method are not visible with the other (derived from Table 5). This, because TSX detects smaller avalanches, while the multiorbital methods detects also avalanches which are otherwise in slopes close to layover.
Finally, the automatic methods detect larger avalanches fairly comparably to the manual method (derived from Table 6), however, weakly visible avalanches and small avalanches which have not been automatically detected cause a miss rate of 435 about 30 %. The apparently lower differentiation of avalanche by manual analysis results from the fact that the automatic method often detects multiple patches instead of a single avalanche (Fig. 9).

Conclusions
We studied the capabilities of the radar satellites TerraSAR-X (TSX) and Sentinel-1 (S1) to detect avalanches in two-image change detection images and multiorbital change detection composites. Manual avalanche mapping results from the high-and 440 medium resolution radar data (TSX, S1) and high resolution optical data (SPOT-6) were compared to each other. An automatic detection method was developed and compared to the manual mapping results.
We conclude that both, TSX and S1 radar images can provide valuable, weather-independent information about avalanche activity, even in difficult alpine terrain. Despite of different lower cut-off sizes of about 500 m 2 for TSX and 2000 m 2 for S1, avalanche outlines and the total area of mapped avalanches agree very well with each other.

445
Between the manual TSX and SPOT-6 mapping results, we found a fair agreement but the total mapped avalanche area of TSX covers only one third (the deposition zone) of the total mapped area (release, path, deposit) in SPOT-6 images. Interestingly, many avalanches located in the cast shadow of SPOT-6 image were not detected whereas they were clearly visible in a TSX image acquired 10 days later. With the automated detection algorithm we found about 60-80% of the avalanches manually mapped in the same image, at least when no large number of old avalanches were present. 450 We found that the non-systematic acquisition program and the possibly high cost can be considered as drawback of TSX data.
Also, with the maximal swath width of 30 km in stripmap mode and a nominal revisiting period of 11 days, an operational use for avalanche mapping over Switzerland with TSX is not feasible. However, the high resolution images can provide valuable information for validation of lower resolution mapping results for pre-defined test sites and if acquisitions are scheduled in advance.
Despite of the lower resolution, we found that S1 provides a convincing solution for systematic avalanche mapping because of the total swath width of 250 km and the revisit period of 6 days when images from the same orbit of both satellites (S1-A and S1-B) are combined. The results from Norway by Eckerstorfer et al. (2018Eckerstorfer et al. ( , 2019 confirm this conclusion.
With the multiorbital combination of systematically available S1 acquisitions from different orbits and with different polarizations we minimized not only areas located in radar layover but also enhanced the radiometric accuracy and obtained a 460 high spatial resolution of about 20 × 20 m. In the resulting change detection image covering entire Switzerland we manually counted in total 7361 new avalanches which occurred during an extreme avalanche period around January 4 th 2018. However, we suppose that mainly avalanches reaching below the wet snow line were detected and that likely many dry snow avalanches were missed because of their lower contrast to the surrounding snow. A disadvantage of the multiorbital composite is the loss of precise timing of avalanches. For operational applications we suggest therefore to minimize the ratio of elapsed time between 465 ascending and descending acquisitions and of the revisit time.
We think that avalanche mapping can be further improved with more advanced methods to combine different orbits, for example with local resolution weighing, LRW (Small, 2012). With that, slopes facing off the radar are weighted stronger which does not only enhance the resolution but should also increase the avalanche visibility: we think that the more omnidirectional scattering of the rough avalanche surface dominates the scattering of smooth snow only for slopes facing off the radar. We 470 found that avalanches are hardly visible in slopes facing the radar (but not yet in layover). As with LRW mountain slopes are unequally weighted the probability that avalanches occur between two averaged images is reduced.
Although we could show that radar change detection mapping with TSX provides results comparable to optical SPOT-6 direct mapping, we note that our study focuses on the exceptionally warm January 2018 with frequent surface melt but also with very intense snowfall periods. As the relative brightness of avalanches with respect to the surrounding snow depends on 475 the water content and the amount of deposited snow, avalanches might be less visible during cold weather with little snowfall. Therefore, we think that an analysis of longer time series of radar based avalanche mapping will provide insight how snow and weather conditions affect the detection rate of radar based methods.
Data availability. TerraSAR-X data are available from the archive https://terrasar-x-archive.terrasar.com. Copernicus Sentinel-1 data processed by ESA have been downloaded from the Copernicus Open Access Hub: https://scihub.copernicus.eu and from the Alaska SAR Facility ASF DAAC 2018 https://www.asf.alaska.edu. The manual mapping results from the optical data and the Sentinel-1 change detection composite of Switzerland is available online Leinss et al., 2019). Table A1. List of S1 acquisitions used for the multiorbital RGB change detection composite shown in Fig. 7 Figure A1. SNAP workflow to process S1 data. The red dashed box is used for creation of the layover and shadow map.  SLF: Wochenbericht 26. Januar-01. Februar 2018, https://www.slf.ch/de/lawinenbulletin-und-schneesituation/wochen-und-winterberichte/