The measurements of snow cover volume in the area of Szrenicki Cirque

The purpose of this project is to determine a fast way of calculating the volume and distribution of snow mantle, which is located in wide terrain concavities in mountain areas. Our study area was so-called Szrenicki Cirque (Kocioł Szrenicki), which is the nival recess, located in Karkonosze Mountains, Poland. We analyzed modern technologies, that are designed to generate 3D-models: terrain laser scanning and close-range photogrammetry (including structure-from-motion technique). There were two major reasons for our research. First, analyzing if a structure-frommotion based software is capable of creating a 3D model of snow cover since potential tie points for adjoin pictures are scarce. The second was to establish the quality and accuracy of this model in relation to potentially more accurate terrestrial laser scanning results. An important issue was also is to estimate the fastest, simplest and least expensive methodology that can be implemented as a daily task of Karkonosze National Park workers. A proper, fast, safe and accurate method of calculating the snow cover volume would incise the safety and avalanche risk evaluation in the vicinity of Karkonosze Mountains. In addition, the developed method can be used to monitor the risk of local spring floods.


Introduction
Estimation of water resources stored within snow cover is important for climatological, water resource management and meteorological applications.Taking into consideration rapid climate changes being observed nowadays it has become apparent that detailed monitoring of snow cover changes might be crucial for hydrological safety of mountain areas [1].The current approach estimates the amount of water contained by the snow cover with Snow Water Equivalent (SWE), the snow cover thickness or by its volume.Those indicators are usually calculated in global or regional scale.However, for local scale approaches close range surveys and measurements are being performed as well.Currently they are being done with the use of GPS, interferometric reflectometers etc. [2], [3].Also due to rapid changes in all survey methods new approaches are emerging namely Terrestrial Laser Scanning.Laser scanners need to be customized due to the fact that snow and ice has very specific reflective abilities.Scans are used either as standalone methods as complementary methods for larger scale analysis done with the use of satellite base remote sensing methods, airborne Lidar and other platforms using Imaging Spectrometer (Airborne Snow Observatory) [2]- [8].Also photogrammetric methods, are being used in order to complement larger scale data [9].This includes standard stereo pictures and newer structure-from-motion methods.
The aim of this study was to check the possible use of new photogrammetric methods, structure-from-motion algorithms with pictures done with the use of hand held noncalibrated digital camera in estimating the thickness of snow cover.The scope of this study was to check the accuracy of photogrammetric model and then calculate the volume of the snow cover as a difference between models with or without snow.This will lead to providing the information on snow cover thickness in every point of the study area.The total difference could lead to estimating Snow Water Equivalent.

Study area
Study area is so-called Szrenicki Cirque (Kocioł Szrenicki) a nival recess, located in Karkonosze Mountains, Szrenicki Massif (1362 m a.s.l.) Poland.It is placed between Szrenica mountain and Szkolnik mountain (fig.1.).From administrative point of view it is located in Lower Silesia Province, in Jelenia Góra County, in Szklarska Poręba gmin [10], [11].From geomorphological point of view this object is a nival recess, however it is still called a "cirque" (Kocioł Szrenicki), due to the fact that it has similar appearance to real Karkonosze Cirques located in close vicinity.Karkonosze are the highest mountain chain of Sudety mountains.They are located at the Polish -Czech Republic border.They have about 30 km in length and 20 km in width.Its highest pick is Śnieżka mountain (1602 m.a.s.l).Karkonosze are some of the youngest mountains that were mostly formed during Alpine orogeny phase along with all of the Bohemian Massif.Its unique morphology is due to double orogeny phases that took part in forming process.Sudety were first formed during Hercynian orogeny and then eroded to the point when no visible mountain forms occurred.Later Alpine orogeny formed mountains that are observed today but since it had to go through leftover parts of the first forming process it resulted in unique Fault Block Moutains that can be observed in some parts of Sudety [10] and are not observed in the Carpathian or Alps Mountains.
Snow cover in the vicinity of Karkonosze mountains is abundant and long lasting (it stays much longer that in other similar regions) [11].Multiple researches on its nature were published thus far [12].The cover size and amount of snow material varies largely with local climate and geometrical conditions.Larger snow covers were observed in the area of cirque and a nival recesses.This includes also Szrenicki cirque, where snow cover can be found for as long as May/June.There is little to no evidence of permafrost occurring currently in Karkonosze [13].
Study area, Szrenicki Cirque was chosen due to the fact that the changes in snow cover, leading to water table seasonal upraises, are being currently observed by the authorities of the Karkonosze National Park.Also its closeness to tourist roads and its uncovered structure was a factor.The study was performed in cooperation with The Base Stations of the Integrated Environmental Monitoring KARKONOSZE [14].This institution also performs survey on physical and chemical properties of snow cover.

Survey methodology
In order to gather data it was necessary to perform two survey seasons.They were done with the use of two current survey methods.First method included using terrestrial laser scanning.It was done with Leica C10 [15].This method provides a point cloud, discreet spatial information on scanned objects.Point cloud consists of a large number of points measured on the surface of the object.Besides XYZ coordinates a point cloud caries information on intensity (strength of the laser signal coming back to the instrument) and RGB colour at the measured point.In order to perform registration, connecting two adjoin scans, it is needed to define the same points on them.This can be done with natural, easy to identify elements (window edges, plat surfaces, nails etc.).However this approach does not work with natural, topographic surfaces since they do not possess such features.This is why artificial points -targets, need to be placed on the object.Targets are usually black and white chessboards, ideal spheres or, in case of Leica scanners HDS (round blue) elements.In order to place point cloud within global coordinates, targets XYZ coordinates need to be measured and results inputted during registration process.
Second survey method includes new photogrammetric algorithm called structure-frommotion. 3D model, in the shape of point cloud or grid models, is created from a group of images done by hand held non-calibrated camera [16], [17].Structure-from-motion algorithms allow the recognition of high-contrast features (of the object), following their movement through series of pictures and produce a sparse point cloud basing on features placement in the image series [16].In most current instalments cameras do not need to be pre-calibrated, camera calibration is done within the algorithm.Created models are accurate in regard to the geometric features but have an incorrect scale.In order to add scale and global coordinates, similarly as with laser scanners, easy to identify points of known coordinates must be given.In case of photogrammetric survey those points are almost always black and white chessboards [16], [17].
In order to calculate snow cover volume two survey seasons were performed.First, in the early spring of 2017 when snow cover was at maximum a photogrammetric survey season was performed.It was done while moving along a baseline located nearly parallel to the object.Complementary parts of the object were visible on at least 3 images, that covered the object top to bottom.They were taken while aiming optical axis parallel to each other.The procedure was repeated with other baselines (also nearly parallel to the object) but pictures on each line were taken from differed angles (of optical axis).During each season axes needed to be parallel in horizontal and vertical aspect.There were 5 targets placed on the snow cover and there coordinates measured with the use of GPS receiver.Targets were placed on differed heights and no 3 of them created a line.
The second survey season was done in the early autumn of 2017 it consisted of two measurements done simultaneously.Laser scanning survey consisted of 4 scan stations that were registered with the use of HDS Leica targets, chessboards and spheres.Survey season was taken during extremely unfavourable conditions, during heavy rain but results were satisfying.Also a photogrammetric season was done in a similar way as in winter.In order to register scans, scale picture based model and provide global coordinate system a total of 12 targets were placed on the object, 5 of them were measured with the use of GPS receivers.The second season provided with not only data for digital terrain model but also was used to estimate the accuracy of structure-from-motion algorithms.Since those algorithms are based on automatic edge search, and topographic objects do not provide many of those performance in those conditions needed to be checked.

Data pre-processing
Point clouds obtained during the laser scanning season were connected with the use of Leica Cyclone software.Registration results were at 7 mm of total accuracy.This might be viewed as low accuracy but since the field conditions were so difficult and displacement of targets possible those results were regarded as good.
Pictures were calculated with Agisoft Photoscan software that is known to be based on traditional stereo-image analysis and performs best in objects exterior modelling.Agisoft automatically sorts uploaded pictures into smaller groups that consist of images done with the same camera parameters, focal length etc. [17] The next process, Align photos, allows finding the same elements on two or more images.This results in a sparse point cloud of complementary points and visualisation of estimated camera location.Next steps Build dense point cloud and Build mesh allow for interpolation of existing data into a coloured point cloud and colour blended mesh.A 3D geometric model is being created, with incorrect scale.In order to create a 3D metric model or a 3D georeferenced scaled model it is necessary to add the targets' coordinates.The result accuracy calculated by the software was at 10 cm.As a result 3 point clouds in the same coordinate system were created.

Volume calculations and 3D models analysis
After creating comparable point clouds of appropriate precision and the same coordinate system it was possible to perform 3D analysis.Next, the most time taking step, due to needed computer capabilities, is creating terrain models.In order for MESH-generating algorithms to perform correctly it is necessary to clean point clouds from all points not representing terrain; this includes artifacts, rocks, people, trees etc.The cleaning process was done in Bentley Pointools software, which turned out to be very convenient due to its low hardware requirements and multiple navigation options.Groups of points representing trees, shrubs, artificial elements in the outskirts of point clouds were removed.One can only create a correct Digital Terrain Models (DTM) based on points representing topographic surface.
After preparation three DTM's were built.For the purpose of this paper, the following identifications were proposed: • "zima_fotka" -a 3D model basing on photogrammetric measurements, which took place in winter, when a snow cover still lingered (due to some problems with reading big data in programs QGIS and Cloud Compare, it was necessary to subsample this point cloud from about 36 million to about 12 million points and it was done in Cloud Compare) (In QGIS and Cloud Compare subsampled version was analyzed, in ArcGIS subsampled and full version were analyzed), • "jesien_fotka" -from photogrammetric measurements, which took place in autumn, when there was no snow cover yet, • "jesien_skaning" -from Terrestrial Laser Scanning (TLS) measurements, which took place also in autumn.
For the purpose of creating DTM's and to performing volume calculations based on them, following applications were used: • ArcGIS Desktop 10.3.1, which enables to create DTM's in the form of TIN and GRID, • QGIS 2.8.13 (Las Palmas), which enables to create DTM's in the form of GRID, • Cloud Compare v2.6.3, which enables to create DTM's in the form of GRID only for (and during) volume calculations.
In order to compare volume calculations and there accuracy resulting from different methods, DTM's have to be built over one, common area on a horizontal plane.Performed surveys covered different areas with usable DTM models.So it was necessary to determine the common area for analyzed DTM's, after creating them with all methods, and then "cut" all these DTM's up to this common area.At this point it is worth mentioning that the other possible way of cleaning point clouds from points located out of the common area exists.DTM's are then generated based on these point clouds.This solution would be faster, because the time needed to generate DTM's primarily depends on size of input data, but it has its disadvantages: • it's very difficult to interpret manually a common area, covered by point clouds, in programs used for editing point clouds, • while generating DTM's with interpolation methods for areas near the boundary of the common area, it is appropriate to use also points located out of this boundary, also representing terrain correctly.This way the model is more accurate within.Therefore, this solution was never an option, with the exception of volume calculations in Cloud Compare.In this case it was necessary to use point clouds after "cutting" them up to common area of DTM's used in ArcGIS and QGIS.Editing DTM's in Cloud Compare is not possible.
DTM's created in ArcGIS and QGIS were cut up to common areas with the use of two options: • common area for DTM's based on all point clouds, i. e. "zima_fotka", "jesien_fotka", "zima_skaning" (purpose: estimating compliance of both survey methods and of all used computer programs), • common area for DTM's based on point clouds from photogrammetric measurements, i. e.
"zima_fotka", "jesien_fotka" (purpose: estimating compliance of computer programs with each other, for research material from one method).
In order to perform a 3D analysis two programs were used.Their descriptions are given below.ArcGIS software reads point clouds in LAS format, so first it was necessary to export them to this format (it was done in Cloud Compare).In ArcGIS the specific sets of data, called "LAS Datasets", readable only in this program, were created from the input point clouds.This format allowed the creation of two forms of DTM.First TIN's with setting: random point thinning with the maximum number of nodes set at 2,000,000, then GRID's with linear interpolation (so-called triangulation) method for estimating elevation of nodes, with no point thinning, with mesh size equal to 0,25 m was created.DTM's designed this way are readable only for ArcGIS.QGIS software on the other hand reads point clouds in text formats, e. g.PTS.Only co-called GRID's were created.Estimating elevation of nodes was done with triangulation method, mesh size was at 0.25 m, similar as in ArcGIS, with only one exception: GRID based on the point cloud "zima_fotka" had mesh size equal to 0.50 m.
Volume calculations in ArcGIS and QGIS were realized only after determining the common area and after reducing extent of DTM's to these areas.The calculations consists of assuming some horizontal plane with established elevation, definitely under whole DTM's, and then volumes between model and this plane were calculated.To determine volume between two models, appropriate results were subtracted.Results of calculations realized with these programs are presented in tab.1.: In the program Cloud Compare volume calculations were realized in two ways: • checking differences between two point clouds (basing on GRID DTM'S -mesh sizes: 2 m, 1 m, 0.25 m), • checking differences between theoretical "base plane" defined under all other point clouds and DTMs -calculation for mesh size equal to 0.25 m. Results of calculations realized with this program are also presented in tab. 1.
There is a possibility to make series of analysis in this application basing on imported point clouds.To illustrate spatial characteristics of analyzed point clouds, the following analysis were realized: • examination of vertical distances between following point clouds: "jesien_fotka" to "zima_fotka" and "jesien_skaning" to "zima_fotka" (for recognition snow cover thickness distribution) and "jesien_skaning" to "jesien_fotka" (for comparison of difference between point clouds from both autumn survey methods) -the cloud "zima_fotka" in subsampled version; • examination of points density in all three point clouds (in each point of cloud points within 0.12 m radius were counted and divided by area of this circle, results in pts/m 2 )point cloud "zima_fotka" in subsampled version.

Discussion and results
The table below presents all the results of previously described analysis.It allows for better interpretation of volume results from the softwares used.The table 1 presents the volume calculations.In bold one can see the maximum number, which was given 100% arbitrarily.All other calculations were done in relation to this maximum.When taking into consideration differences between spring -autumn models all results seems to be similar.In regard to maximum value those results do not exceed 10% percent.This result occurs when the analysed region is limited to common area.When the common area is limited to two models (zima_fotka, jesien_fotka) differences are nearing 20%.This might be a result of errors of structure-from-motion algorithms occurring on the area edges, errors in mesh generating procedures and others.
Surprisingly differed results occur in data prepared and processed in Cloud Compare software.Analyses of results from one surveys season reach up to 70%.Maximum is at average 1 m GRID size and minimum at the largest 2 m GRID size.Those errors might result from Cloud Compare algorithms programming that might be not suitable for such objects.Even though percentage differences are big, it is worth mentioning, that absolute values of these differences on such wide area are relatively small.The above mentioned table also shows that decrease in the point count in the zima_fotka point cloud does not indicate significant negative volume calculation results.Table 2 presents volume differences between jesien_zima and jesien_fotka or jesien_skaniang.All differences were calculated for the same limited area that is covered by all models.The reference difference was zima_fotka -jesien_skaning.Percentage differences analysis shows that in case of ArcGIS those numbers are similarly small and snow volume is only a little larger if jesien_skaning was used as a reference model, in case of QGIS and Cloud Compare those numbers are also small, but if jesien_skaning was used as a reference model, snow volume is a little smaller.
The summary above describes a series of experiments to show that differences between 3D models obtained from laser scanning and structure-from-motion algorithms give similar results.This implies that those approaches can both be used for modeling topographic surfaces, and since the photogrammetric method is much faster and requires less financial input it is highly favorable in this case.Also differences between commercial ArcGIS software and free QGIS and Cloud Compare software are relatively small, so from a practical point of view it allows the reduction of costs of such calculations.Those preliminary results show a lot of promise in using surface based structure-from-motion algorithms in calculating the volume of snow cover.However further investigation into analyzing accuracy of winter models are necessary.
Analysis presented in this paper may be important not only for scientific reasons, but also are response to the needs related to tasks of Karkonoski National Park employees.Presented methodology can be used for cheap, fast measurements of snow cover volume.

Fig. 1 .
Fig. 1.Study area a) Karkonosze placed on polish administrative province maps, b) the Karkonosze mesoregion in the background of the division of Kondracki physiogeographic provinces (the Bohemian Massif is marked gray), c) approximate location of the study area on the background of the Karkonosze mesoregion and the Polish border, d) study area on the topographic map 1: 25,000 (Main Land Surveyor).The range of the map visible on part d is marked with a rectangle in part c.

Table 1 .
Summary of the results of volume calculations [m 3 ].

Table 2 .
Comparison of snow cover volume calculations results using two survey methods in the autumn season [m 3 ].