UAV-BASED LIDAR HIGH-RESOLUTION SNOW DEPTH MAPPING IN THE SWISS ALPS: COMPARING FLAT AND STEEP FORESTS

Snow depth mapping in Alpine forests is of high importance for hydrogeology, ecology, tourism, and natural hazards prevention. Different remote sensing approaches have been employed for the precise mapping of snow depth within forests. However, optical sensors cannot provide below-canopy information. While Airborne Laser Scanning (ALS) systems have been used successfully in this context and allow obtaining data below canopies, the costs of acquisitions are very high, not allowing frequent data acquisitions. UAVbased Lidar technology potentially can provide the critical below-canopy information at lower cost and allows for frequent acquisitions. First attempts to employ a UAV-based Lidar system in forests have proven promising, but they are limited to flat forests and to gridlevel snow depth calculations. In this study, we present UAV-based Lidar data of both flat and steep forests. Different Lidar processing workflows are analyzed and compared, and snow depth algorithms are used both at the point and the grid level. Whereas the UAVLidar system proved capable of mapping snow in both environments, the steep forests' data processing comes with greater challenges, especially for the 3D registration, ground classification, and point-to-point snow depth calculations.


INTRODUCTION
Forests cover around 30% of the surface of Switzerland. Especially in Alpine regions, the snow dynamics and the snowpack properties affect ecology (Wipf et al., 2009), water resources (Farinotti et al., 2012;Thornton et al., 2021) and are relevant for snowmelt forecasts and tourism (Pütz et al., 2011). They are also important in the context of natural hazards such as avalanche prevention (Einhorn et al., 2015).
However, mapping snow in forests remains a big challenge for remote sensing applications. Unmanned Aerial Vehicle (UAV)mounted optical sensors appear inefficient in mapping snow below canopies (Harder et al., 2020). Airborne Lidar Systems (ALS) show promising capabilities in mapping snow in forests within different scales and canopy types (Hopkinson et al., 2004;Currier et al., 2019). However, apart from the high cost of data acquisition, challenges related to the forested environment can underestimate the snow depths (Hopkinson et al., 2004). Also, data gaps below the canopy (Mazzotti et al., 2019;Broxton et al., 2015) can limit ALS systems' application to map snow dynamics. UAV-based Lidar techniques potentially overcome the constraints of optical sensors by capturing both canopies and the sub-canopy snow surface at the same time. The feasibility of frequent scans of various scanning angles in contrast to the ALS systems leads to denser point clouds and differential snow depth grids of high spatial and temporal resolution. Only a few UAVbased Lidar studies for forest snow mapping have been published (Harder et al., 2020;Jacobs et al., 2021). The challenges identified include accurate flight planning and the subsequent Lidar-related processing workflows. These challenges are accentuated in steep forested terrain due to the difficulties of operating UAV in such terrain. However, steep terrain typically leads to increased vertical inaccuracies for Lidar systems with larger laser beam footprints (Deems et al., 2013;Baltsavias et al., 1999). For this reason, UAV-based Lidar has the potential to * Corresponding author generate snow cover datasets of higher quality than those that can be obtained by ALS in steep forested terrain, yet such an application has not been attempted yet.
In this study, we present a UAV-based Lidar approach and compare snow mapping techniques over flat and steep forests. For the first time, we present UAV-based Lidar data obtained over steep slopes with a high resolution of 0.5m and forest snow depth extraction directly from the point clouds. The study's specific objectives are 1) to present Lidar data processing workflows for snow cover mapping in flat and steep forested terrain, and 2) to compare 3D registration algorithms and the snow depth mapping techniques applied individually to each case study.

UAV acquisitions
UAV-based Lidar data were acquired over a flat forest of 200 by 200 m at Flin, Engadine, in the Eastern Swiss Alps in 2019. One snow-on data acquisition in late March and a snow-off acquisition in the middle of June were carried out. In 2020, data over two forested slopes of opposing orientation were collected over an area of the same extent in Flüela valley near Davos. Eight and thirteen snow-on flights over the north and south-exposed slopes were complemented by snow-off flights at the end of April and at the beginning of May, respectively.
We used a Yellowscan Mapper II mounted on multicopter from Altigator SA. The Lidar system operates in the 905nm wavelength with a FOV of 100• and horizontal and vertical accuracies of 15 and 5 cm, respectively, and weighs approximately 3 kg (including battery). The multicopter can carry a payload of 12 kg.
The flight planning for both 2019 and 2020 combined four interlaced flying patterns for each acquisition to achieve higher point densities. Each flying pattern consisted of four scanning lines. The flying altitude was retained at 70-80 m above ground, aiming to follow the topography in the steep forest case, and the flight speed was set to 3m/s. These flying characteristics led to a mean point density of 30 points/m 2 per scan with a maximum of two returns per pulse.

In-situ snow depth measurements
In-situ snow depth measurements were acquired for both studies to validate the snow depth maps. In the flat forest case, snow probe measurements were acquired in five different transects within the forest with a one-meter spacing. The measurements took place after the snow-on flight so as not to disturb the snow surface. 74 snow probe measurements were acquired in total within the five transects.
In the case of the steep forests, due to the difficulty in the accessibility of the sites and the need for in-situ snow measurements for each winter acquisition, we followed different procedures: We used five snow depth stakes that we installed at fixed locations within the forests. They remained in the fields during the whole winter and were removed after snowmelt. In addition to the snow depth measurements on the snow depth stakes, we placed three plastic hemispheres covered with reflective tape on the snow surface during every acquisition. The reflective tape could backscatter the incident beams of the Lidar in the wavelength of 905 nm. Therefore, the hemispheres could be identified on the point clouds (PCs) using the "intensity" value of each point that is automatically recorded by the Lidar system and saved as a scalar field on the PC. Absolute snow depths with a snow probe were measured in the centre of the hemispheres. The snow depth stakes and the hemispheres used provided 8 insitu snow measurements per acquisition.

Post Processing Kinematic correction of the UAV trajectory
The data from the inertial navigation system of the Lidar system that include both the IMU and GPS recordings were retrieved and converted into PCs using Yellowscan CloudStation (version 1.18.0, licensed). For increased PC accuracy, the IMU and GPS data were further used to correct the raw UAV trajectory. This integrated forward and backward GNSS post-processing corrections, which were performed using the POSPac UAV software (version 8.3, licensed). We used the GNSS data from a base station from the Swiss Continuously Operating Reference Stations (CORS) network within 2km of the field sites. The Post Processing Kinematic (PPK) correction of the UAV trajectory led to position accuracies around 1-2 cm. The corrected trajectory was then exported as a Smoothed Best Estimated Trajectory (SBET) file and re-imported into the CloudStation to correct the PCs. This workflow was followed for all the data acquisitions of the flat and steep forests.

Manual strip adjustment
As a next step, for each acquisition, a strip adjustment quality check took place to evaluate the alignment of the four scanning lines. In some cases, even after the SBET-based post-processing of the PCs, there was a misalignment between neighboring scanning lines (hereafter strips). This can be related to uncertainties in calibration angles, unsuccessful initialization of the IMU, satellite loss, or heading drifts. We thus established a manual strip adjustment workflow that was applied in case of misalignments after the SBET-based postprocessing of the PCs. This workflow includes the selection of a random intersection of the two neighboring strips within the PC. The workflow took place on CloudCompare software (version 2.11, open-source). One intersection remained fixed (reference), while the other was manually moved to register to the reference. The 3D transformation matrix was then extracted and applied to the whole strip. A quality check was done to evaluate if the transformation matrix was well applied to all the corresponding strips' points. The process was repeated using different intersections until the extracted transformation matrix was successfully applied to all the points. The workflow was applied separately for every two strips, and the result was one PC for each scan.
The manual strip adjustment was mainly used for the steep forests data, where the higher alignment errors were identified within neighboring strips. For the flat forest, the SBET post-processed PCs yielded better strip alignment. Therefore, it was not needed to apply the manual workflow.

Registration of the scans
The strip adjustment quality check mentioned in section 2.4 was followed by the registration of each acquisition's four flying patterns (hereafter scans). The aim was to extract one PC per acquisition, combining all four scans for increased point density.
For the data of the flat forest, the Iterative Closest Point (ICP) algorithm was used for the registration of the four scans as follows.
For every two scans, the ground points of the overlapping parts were used for the ICP algorithm. The scans' registration per acquisition for the flat forest took place after the snow/ground classification (see section 2.6). The ICP algorithm was applied only to the snow/ground points. The vegetation points were excluded since the canopy point density was not satisfying enough for the ICP algorithm to operate adequately. After applying the ICP algorithm on the ground points of the overlapping parts of the two scans, the 3D transformation matrix was extracted and then applied to the entire PCs of the corresponding scans. The quality check included selection of random intersections within the PCs to evaluate whether the 3D transformation matrix was well applied to all the points. If not, the ICP algorithm was re-applied but in a smaller sub-sample of the points within the scans' overlapping part until the 3D transformation matrix led to satisfactory matching of all the points throughout the whole PCs.
This process took place for every two overlapped scans until all four of them were correctly registered to each other. The final result was one PC per acquisition. The ICP implementation, as well as the quality check, were performed using the CloudCompare software. Interestingly, the ICP algorithm could not be applied to the forested slopes' data as it led to large errors (see section 3.2.2). Alternatively, we performed a manual strip adjustment as described in section 2.4 using the whole scans instead of individual strips.
The registered PCs derived from the application of the ICP and the manual adjustment workflow were further used for the snow/ground classification.

Snow/Ground classification
The registered PCs derived from the application of the ICP and the manual adjustment workflow were further used for the snow/ground classification, which is an essential step for creating the snow depth maps. It aims to separate the snow or ground from the vegetation points in order to extract the snow or bare earth surface. For this step, the lasground tool from LAStools software (academic version 198012, licensed) was used, which classifies points into ground and non-ground points.

Ground classification at the flat forest
In the flat forest case, for the snow-on PC, the lasground tool was used with the default parameters proposed for a flat natural environment (parameter set called "natural"). Additionally, we applied the setting "offset = 0.6 m", which entails all points within this threshold from the initial ground-classified points to be classified as ground as well.
For the snow-off PC classification, the lasground algorithm had to account also for the underlying low vegetation in the middle of June when the snow-off acquisition took place. The lowvegetation points could be wrongly classified as ground points and therefore lead to an over-estimation of ground elevations in the DEM of the snow-off acquisition, which would result in an underestimation of the snow depths. Despite the "natural" and "step 5" parameters suggested for the flat forest, the offset parameter was set to 0.8 m to include all ground points.

Ground classification at the steep forests
The snow/ground classification workflow followed was different from the one used for the flat forest. In the case of the steep forests, the algorithm had to identify whether the differences between points occurred because the points are of a different class (i.e. capturing different surfaces), or due to the elevation difference.
Since the south-exposed site was steeper than the north-exposed one, the vertical inaccuracies derived from the in-built Lidar configuration were higher. This constraint affected the selection of the lasground parameters for the snow-on PCs, which were then different from those used for the north-exposed slope.
Regarding the snow-off PCs' ground classification, the algorithm had to account also for the low-vegetation points. Since the snowoff acquisition for the north-exposed slope took place later than the south-exposed slope, the ground classification parameters used were different.
For the snow classification, for the south-exposed slope, the offset parameter was set to 0.7 m. To remove all the spikes 5 m above or below the initial snow classified surface, we set the "spike" parameter equal to 5 m. To account for the forested environment, we used the parameter "nature," and in order to intensify the search of initial ground points due to the terrain's steepness, we added the parameter "hyper fine". With these selected parameters resulting from many iterations, we managed to separate the snow from the vegetation points even in the forest's steeper parts.
For the snow classification of the north-exposed slope, the parameters remained the same as those for the south-exposed slope. The only difference was the offset that was set to 0.5 instead of 0.7 m.
Regarding the ground classification of the snow-off PCs, the lasground tool had to account for the low-vegetation points and the steepness of the terrain. For the south-exposed slope, we used the same parameters for the snow-on PCs. The offset in the case of the snow-off PC was set to 0.5 m.
In order to remove the low-vegetation points, we additionally applied the lasthin tool. The ground classified points from the lasground tool were used on lasthin to keep only the percentile of 30% of the elevations within a step of 10 cm, similar to Mazzotti et al. (2019). This "percentile 30" parameter aimed at retaining as ground points only the lowest points within the 30 % of the elevation values of all points within a step of 10cm. This concluded to a thinner ground surface that most of the lowvegetation points had been filtered out.
For the north-exposed slope, the lasground parameters were the same as those for the snow classification, with the difference of setting the offset to 0.8 m. As in the south-exposed slope, the "percentile 30" option on lasthin tool was used to filter out the low-vegetation points within a step of 10cm.
To evaluate the classification results for both flat and steep forests, we were visualizing the two classes via the "classification" scalar field assigned to the PCs during the classification process from lasground. Afterward, we took intersection within the classified PCs to evaluate whether the classification was efficient, especially in challenged parts, such as below-canopies. Not satisfying classification, including misclassified vegetation points or fluffy snow surface led us to change the lasground parameters. The evaluation process was repeated until we extracted satisfying classification results, including as less as possible misclassified points and thin summer ground surface.

Snow mapping over the flat forested terrain
The snow and ground classified PC were used for the creation of the snow depth maps. The workflow followed for the flat and steep forest was different.

Cloud to Cloud (C2C) distance algorithm
The distances between the snow-on and snow-off ground classified points were estimated using the Cloud to Cloud (C2C) distance algorithm on CloudCompare. The C2C algorithm uses the one PC as a "reference" and computes the Euclidean distances to the second one that acts as the "compared" PC. The identification of the closest point on the compared PC is based on the Hausdorff algorithm (Rockafellar and Wets, 2005).
Before applying the C2C distance algorithm, we tested the relative position between the snow-on and snow-off PCs. They are dependent on the geolocation inaccuracies between the acquisitions. These are relative to the initial processing (SBET post-processing, strip adjustment, and registration of the scans).
Despite the precise application of the registration method at the strip or scan level (ICP or manual strip adjustment), the registration comes with errors from millimeters to centimeters, depending on the quality of the PC (density, noise) and the accuracy of the method used. Since there are inaccuracies at each processing level, the final PCs retain all the processing steps' errors. Consequently, it is important to assess the geolocation accuracies between the PCs.
We used the groomed street on the field site to detect potential misalignments between snow-on and snow-off PCs. Since the groomed road was snow-free during the snow-on and snow-off acquisitions, the two PCs were expected to be finely registered along this road, which was indeed the case. Consequently, no additional registration was required.
For the C2C distance algorithm, we used the ground classified snow-off PC as reference, and the Euclidean distances to the snow classified snow-on PC were computed. The distances were stored as a new scalar field into the ground classified snow-off PC and gridded into a 0.5m spatial resolution raster to yield a snow depth map. For the gridding, we selected the average values of points falling within a pixel of 0.5 by 0.5m. The gridding was also performed using Cloud Compare without using any type of interpolation for gap filling purposes.

DEM of difference (DoD)
The second methodology, called DEM of difference (DoD) is based on the subtraction of two DEMs to create a new one. The two initial DEMs used for the subtraction are the snow-on and snow-off DEMs. The produced DoD is the snow depth map. This second approach acquires as first step the gridding of the ground points of the snow-on and snow-off PCs. The spatial resolution chosen for our study was 0.5m depending on the point density of the PCs. Finer resolutions of 0.2 and 0.3 m were tested but yielded many data gaps. The PCs were gridded using the lasgrid tool from LAStools. Each pixel was given the average elevation value of all the ground points within an area of 0.5 by 0.5 m. The resulting elevation rasters were the DEMs of the snow-on and snow-off acquisitions. By subtracting the snow-off DEM from the snow-on DEM, we created the DEM of difference (DoD) in which each pixel represents the snow depth value. Figure 3. Snow depth grids extracted from DoD lasgrid and C2C algorithms.

Snow mapping over the steep forested terrain
The snow maps for the steep forested slopes were extracted using the before-mentioned DEM of difference (DoD) method. The C2C algorithm was initially tried on the data but yielded high errors due to the terrain's steepness (see section 3).
The steep forests lacked snow-free features within the sites that would have allowed evaluation of the alignment of the snow-on and snow-off PCs. To assess and minimize geolocation errors, we could only use tree points to register different PCs. As we were flying during wind-free conditions and when the snow intercepted in the tree branches was limited, the trees were stable features within all the acquisitions. The points corresponding to the branches and the treetops were therefore used for applying different feature-based 3D registration algorithms: point-to-point manual alignment, ICP algorithm, and 3D transformation proposed by Ferraz et al. (2018) were applied to our data. The first two point-based methods did not identify pairing points in the trees (more details in section 3). Alternatively, we used the registration workflow proposed by Ferraz et al., (2018). Since this is a registration methodology based on the treetop positions within the scans that need to be registered, we first had to segment the treetops. We used the ITC segmentation package (itcSegment) integrated into R designed by Michele Dalponte (2018) to extract the positions of the treetops on the snow-off and each of the snowon PCs. The 3D optimal rigid transformation was followed based on the extracted rotation and translation matrices from the treetops' coordinates. The snow-off PC was used as reference, and as the registered one was used, the snow-on PC. This method proved to be also not applicable to our data (see section 3.2.2). Taking different transections within the different PCs, it was concluded that the geolocation errors within the different scans The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B3-2021 XXIV ISPRS Congress (2021 edition) existed only in some locations in the PCs' borders. This made the misalignments non-uniform between the PCs. Therefore, the geolocation errors could not be quantified through the transformation matrices that were extracted from the abovementioned registration methods. Considering this datarelated restriction, we proceeded in the snow depth gridding, acknowledging that those errors would be retained in the snow depth maps.

DEM of difference (DoD)
The ground classified points of the snow-on and snow-off PCs were gridded into rasters of 0.5m spatial resolution, using the average elevation value of all the points falling within an area of 0.5 by 0.5 m using the lasgrid tool. The resulting grids were the DEMs of the snow-off acquisition was then subtracted from each of the snow-on DEM and the output DoD raster was the snow depth grid. The same procedure was followed for all the snow-on acquisitions for both slopes. . Snow depth maps for the south (a) and north (b) exposed steep slopes.

Snow depth maps validation
The in-situ measurements acquired as described in section 2.2 were used to validate the extracted snow depths maps.
For the flat forest, a canopy height model (CHM) of 0.5m resolution was applied to the snow-off PC using the lasheight tool from LAStools. The lasheight tool was applied to the snow-off ground classified PC to estimate each point's relative height above ground. Therefore, all the ground points were assigned zero Z coordinate and the vegetation points were assigned as Z coordinate their relative height above ground. This relative height above the ground of the vegetation point represented the canopy height per point and was gridded into a 0.5 m pixel size raster The snow probe measurements' locations could be identified within the CHM by measuring their distances to the canopy edges. To compare point values with pixel values of an area of 0.5 by 0.5 m, we created a buffer of 3 m around the estimated snow probe location/pixel in the CHM. The final snow depth value that we were using to compare with the snow probe measurement at the specific location was the mean value of all the pixels perpendicular to the identified snow probe location/ pixel on the buffer's snow depth maps.
The RMSE between snow probe measured and Lidar-derived snow depth values were calculated for all the 74 locations within the five transects. The RMSE was calculated separately using both the C2C and DoD-derived snow depth maps and the results were compared to evaluate the efficiency of the two algorithms. For the steep forested slopes, a canopy height model was also calculated with 0.5m pixel size using the lasheight tool. The snow depth sticks' position could be identified within the point cloud of the 1 st of April since the hemispheres were placed in the middle of each bamboo stick. Hence, the coordinates of their locations were derived from the PC. The coordinates were afterwards matched with the corresponding pixels on the snow depth maps. Hence, every stick's location was assigned to one pixel within the snow depth maps. The measured value on each stick as well as the measured snow depths on the middle of the hemispheres, were compared with the snow depth value of the corresponding pixel on the snow depth maps. These differences were expressed in terms of RMSE.
Since the steep forests were characterized by more complex accessibility that limited the in-situ measurements as well as by the lack of snow-free areas that could act as a reference as described in section 2.8, we followed an additional validation step. The measured snow depth values on the sticks and on the hemispheres were used to apply a vertical offset to each of the snow depth grids. The offset was calculated as the average difference between in-situ and Lidar-derived snow depths at the corresponding pixels. Figure 5. Schematic workflow of all the processing steps followed until the extraction of the snow depths maps. In blue the algorithms for the flat forest and in green the one for the steep forests.

UAV data acquisitions in forested environments
The UAV operation in the forested Alpine environment proved challenging with respect to the UAV flight planning, piloting and data processing. The flight planning in the steep terrain was much more demanding than that in flat terrain, as it required the adjustment of the flying altitude to the elevation changes at every location. In addition to the challenging UAV operation over steep forested terrain, the flying altitude adjustment forced the UAV to consume more battery capacity, resulting in a smaller surveyed area. In this terrain, restricted site accessibility due to avalanche danger also limited our capability to acquire ground-validation measurements, with further implications for data processing and evaluation (see section 3.3).

Evaluation of the manual strip adjustment
The manual strip adjustment proposed in this study proved to work efficiently for registering overlapping scanning lines and The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B3-2021 XXIV ISPRS Congress (2021 edition) scans from different acquisitions. It can be applied alternatively when registration algorithms such as the ICP algorithm are inefficient due to low point densities, high noise, and limited spatial overlap. Note that this step demands many quality checks in different parts of the registered PCs after applying the 3D transformation matrix. Due to the lower accuracy of the Lidar system in the borders due to smaller overlapping areas and lower point densities, a transformation matrix derived from points in the middle of the PCs may not appropriately represent the borders. Therefore, a quality check including all the parts of the PCs is recommended. The increased error in the strip adjustment on the data steep forests that was observed can be due to the IMU's difficulty operating accurately when the UAV was changing elevations within two neighboring scanning lines to follow the topography.

Evaluation of the PC registration algorithms on the steep forests
As mentioned in section 2.8, both feature-based coarse and fine registration methods were used for the steep forests based on the tree points. Regarding the feature-based methods, the forested environments that were snow-covered and the lack of humanmade constructions did not allow us to apply feature-based algorithms that are using standard features as tie points (e.g., corners of the buildings). Therefore, we were restricted to the use of the tree points for the implementation of any registration algorithm. The point-to-point manual registration based on tree points selected by the user proved unable to assign a uniform transformation matrix to the whole PC.
Regarding the fine registration algorithms, the ICP algorithm applied on the tree points could not effectively perform on our data. The ICP algorithm is highly dependent on the accuracy of the PCs in terms of noise and point densities. The ICP algorithm's inefficiency in operating on our data can be due to the low point density in the trees. This means that the algorithm could not correctly identify the pairing points on the trees due to the beam's very low possibilities to hit the exact same point at every scan. This was also why the manual point-to-point registration method did not work sufficiently well: It was difficult to identify precisely the same points in both PCs.
The rigid 3D registration workflow proposed by Ferraz et al., (2018) was also not applicable to our data, despite its promising application in cases where stable features on the surface are missing. The vectors between corresponding points in the two scans did not have the same orientation through the whole PCs. This caused the inability of the used singular value decomposition (SVD)-based least-square fitting approach to extract a global rotation and translation matrix for the whole PC. The geolocation inaccuracies that we tried to quantify and minimize through scans were not uniform. The non-uniform inaccuracies within scans also restricted us from extracting differential snow depth grids. The mismatches in the relative positions between two consecutive scans were higher in the borders but difficult to quantify and eliminate. The nonuniformity in the inaccuracies between scans can be due to the decreased accuracy of the Lidar IMU in the steep forests.

Snow depths maps validation
The snow-free pixels corresponding to the groomed road allowed the evaluation of the snow depths maps at the flat forest site. However, the road was occupying only a part on the east side of the forest. For the remaining areas, no snow-free features were present. Despite the snow probe measurements in the middle of the forest, there were no snow-free features in the west part of the field site, neither snow probe measurements were acquired. Therefore, the snow depth map quality could not be evaluated within the whole area of the project. That leads to the conclusion that the intensive in-situ measurements and the snow-free areas that can act as a reference for the extraction of the snow depths are of high importance in snow depth mapping studies. Additionally, the soil during June was wet due to the already started snow melting. This introduced a systematic error in the snow probe measurements due to the probe's penetration into the soil. This error could not be quantified and therefore subtracted from the snow depth values as we did in the case of the snow depth stakes. Yet, despite the use of the average snow depth values of all the pixels within the 3m buffer to achieve the best representation of the snow probe measurements, the comparison between values in specific locations with the averaged values of an area of almost 3 m 2 introduced inaccuracies.
The same limitations in the validation process occurred in the steep forested slopes, mainly because of the limited number of validation measurements (only 8 per acquisition) due to the accessibility restrictions. The inaccuracies introduced by the comparison between snow depth values at specific locations with areas of 0.25m 2 were more substantial in the steep than in the flat forest due to the complex micro-topography. Therefore, we placed the snow depth sticks in comparably flat locations with limited micro-topography, so that the snow depths were expected to be uniform around these areas.

Snow/ground classification
While the snow/ground classification was revealed to be one of the major challenges for both of our studies, the classification in the case of the steep forests was significantly more complicated. This was due to the integration of the increased vertical inaccuracies of the Lidar system on top of the IMU performance. In the case of the steep forests, the lasground algorithm had to deal with the elevation changes within each step. More precisely, it had to identify whether the different elevations for points within the same cell are due to the steepness of the terrain or because they belong to different classes (ground/snow or vegetation). The IMU-related inaccuracies that were increased in the steep forests, despite the ground classification procedure, led to fluffy snow and ground surface. To address this, a thinner snow/ground surface is recommended to be created to include only the snow points into a "median" snow surface, as proposed by Martin Isenburg (2021). This workflow is analysed by Koutantou et al. (in preparation), where we present a snow/ground classification workflow for reducing the fluffiness of the snow surface, using the same data over the steep forests. For the flat forest, the classified snow/ground surface was characterized by much less fluffiness. This can be associated with the better IMU performance, already noticed on the strip adjustment processing step. Therefore no extra processing step, including the thinning of the classified ground/snow surface was applied to the data of the flat forest.

The case study of the flat forest
As presented in figure 6, based on the R 2 values, both C2C and DoD methods resulted in snow depths close to the in-situ measured. The RMSE for all the 74 points within the five transects is between 10-19 cm for the C2C and between 5-18 cm for the DoD method. We see that in most transects, the lasgrid snow depths are closer to the in-situ ones. The C2C distance algorithm calculates Euclidean distances between identified pairing neighboring points. The Euclidean distances differ from the actual (vertical) distances that correspond to snow depths, and this explains the differences between in-situ and C2C-derived snow depths. Low densities can cause the identification of wrong neighboring points.

Steep forests
The C2C distance algorithm yielded high errors in the case of the steep forests. The data in the steep forests have increased low noise in addition to low point densities. Additionally, due to the terrain's steepness, the differences between estimated (Euclidean) and actual distances are even higher (figure 8). Alternatively, we used the local modeling technique on CloudCompare on the reference-snow-off PC to have a better model around the identified nearest point and extract distances closer to the actual (vertical) ones. Similarly, we also applied the cloud-to-mesh (C2M) distance algorithm on CloudCompare. Both the local modeling and the mesh-based methods showed high errors especially in the forests' steeper parts and areas with lower densities, such as close to the canopies. More advanced point-based algorithms, such as multiscale model-to-model cloud comparison (M3C2) have not been tested yet in steep forested terrain to our knowledge A first attempt to apply the M3C2 plugin on CloudCompare in our data was made with the default parameters. The algorithm's performance was similar to the performances of the C2C, local modeling and point-to-mesh methods. A more advanced application of the M3C2 algorithm still needs to be tested and evaluated in steep forests.
Therefore, the grid-based snow depth extraction in steep terrain is more robust than point-based algorithms. Depending on the selected gridding resolution and the gridding value (i.e., an average of all elevations), it can lead to satisfying results. The DoD technique is the most common one for snow mapping studies, and our conclusion is compatible with studies already conducted in open sites. Lower resolutions could also eliminate part of the noise during the snow/ground classification process.
Validation measurements of the absolute and relative snow depth values between Lidar-derived snow depth grids and in-situ measured showed an RMSE of 4-16cm for both steep sites (figure 9) before implementing the offset. The implemented offset that was used reduced the RMSE at every snow depth stake location, as expected. Figure 9. Snow depths before and after the application of the offset for the south (a) and north (b) exposed slopes.
The snow depth maps for the two forested slopes provide interesting data to investigate the distribution of the snow patterns relative to the canopy structure and topography (Koutantou et al., in preparation).

CONCLUSIONS
This study presented the first application of UAV-based lidar specifically aimed at assessing its performance for mapping snow depth on flat versus sloped forested terrain. UAV-based Lidar proved to be an efficient technology in capturing snow depths below canopies in flat and steep forests with a fine-scale of 0.5 m. Mapping steep forests included more challenges, from the accessibility of the sites to the data post-processing workflows at the point cloud and grid levels. Higher vertical inaccuracies are identified for the data over steep forests even after the PPK postprocessing of the IMU and GPS data. The lower quality of the IMU performance in the steep forests imposed limitations in the strip adjustment and in the snow/ground classification processes.
Comparing the 3D registration algorithms used, the manual strip adjustment introduced in our study worked efficiently in both terrains. Coarse and fine registration algorithms that were applied alternatively on the treetops of the steep forests were inefficient due to lower point densities and low noise that made the identification of pairing points very difficult. Regarding the snow depth mapping, both C2C and DoD methods worked efficiently in the flat forest, with the DoD performing better. Evaluating the C2C algorithm in steep terrain for the first time yielded high errors due to the above-mentioned data restrictions (low densities and surface fluffiness). Thus, the grid subtraction method (DoD) proved to be more robust over steep terrain. Overall, our results confirmed UAV lidar as a promising tool for forest snow mapping in complex forested terrain, but they also revealed that many technical challenges still need to be overcome.