Application of ALS data for GNSS terrain obstacles inventory

. Current constellation of global navigation satellite system (GNSS) ensures signal availability even in severe observational conditions like urban canyon or under tree canopy. However, positioning in such environment remains a challenge because obstacles can block, reflect and diffract GNSS signals which significantly affects accuracy. Those errors are strongly sight dependent and cannot be mitigated in differential positioning that is why, knowledge of the shape and spatial distribution of terrain obstacles is essential. In this paper using of airborne laser scanning (ALS) data for terrain obstacles inventory is presented and evaluated. In proposed method terrain obstacle models have been derived from ASCII ALS data file using open source QGIS with LAStools software suite and dedicated ALSObstModel plugin. Test models were developed for three geodetic control points with different environmental characteristics. For each point reference model from direct tachometry measurements have been obtained. An average error in determining the elevation of the terrain obstacles from ALS based models was 0.6° to 1.7°. This distance corresponds to 3 to 6 minutes of satellite in orbit.


Introduction
Wide availability of networks of permanent reference stations and current constellation of global navigation satellite system (GNSS), consists of almost 100 satellites make satellite positioning techniques very popular and readily used. The main limitation of these techniques is the need to provide an access to unobstructed horizon. Some terrain obstacles, like in case of urban canyon [1] or tree canopy [2] can block, reflect and diffract GNSS signals, which significantly affects accuracy. Due to the fact that these errors are strongly sight dependent and cannot be mitigated in differential positioning, knowledge of the shape and spatial distribution of terrain obstacles is essential in proper mission planning as well as GNSS data processing. Regardless of whether we are considering GNSS mission planning or pre-processing GNSS data, we always need to relate the position of the satellites to the terrain obstacle layout. Thus, in the procedure, we separate two steps: -creation a polar plot of the elevations of the obstacles above horizon -computation the satellite positions and comparison obstacles' and satellites' elevations [3].
Over the last 20 years, many methods of generating obstacles diagrams have been proposed. These methods can be divided into two categories: first requiring some observations directly at the measuring point, and second one based on spatial data sets. The most obvious method is to use theodolite. This method characterized by the highest accuracy but is very time-consuming and poses some technical difficulties in the case of curtains on high elevations. The more advanced solutions, successfully tested by authors in previous research, are terrain obstacles modelling based on terrestrial laser scanning and digital hemispherical photography [4,5]. The first group of methods allows to obtain very accurate models of terrain obstacles. And it is recommended to use such methods if there is a possibility of quick access to measurement points. In the case of measurement campaigns taking place in large areas, performing direct measurements would consume too much time and resources. Therefore, in such cases, one should look for methods that allows to obtain the appropriate and accurate terrain obstacle models without having to perform a direct measurements.
Therefore, there is a second group of satellite visibility determining methods that includes the most advanced methods based on 3D digital surface models [6] and 3D buildings models. Moreover, 3D buildings models can be used not only for advanced GNSS mission planning, but also for supporting GNSS navigation in urban canyon like in GNSS shadow matching technique [7]. The use of 3D models is, however, burdened with errors in measurement data, model execution errors and the error space model location itself. In addition the main limitation of those methods is the fact, that there are no complete 3D surface models (models including, apart from buildings, trees and other obstacles) for many areas, except large urban agglomerations. Therefore, using 3D models data in forested areas is virtually impossible.
In recent years, LiDAR (light detection and ranging) technology has become a popular method of spatial data obtaining and ALS data become more accessible than ever. For example in POLAND in the years 2011-2014, the LIDAR data for 92% of the country's area was obtained as part of the ISOK project [8]. These data were obtained regardless of the level of urbanization, therefore it is possible to use them for all points of the geodetic network. From the point of view of terrain obstacle modelling, ALS has a number of advantages, especially provides high dense 3D points cloud in which points represents land surface as well as every objects above it.

Object and Source Data
In order to conduct surveys concerning the inventory of terrain obstacles, three GNSS measurement points were selected. The test points were located at the detailed geodetic control network points with catalogue numbers: 0838806, 0839104 and 0838904. For each of the points, catalog coordinates in the flat coordinate system PL-2000 were obtained. All points were located in the campus of the University of Warmia and Mazury in Olsztyn in Kortowo are. When selecting measuring objects, the most important thing was the inventory of various types of terrain obstacles, mainly buildings and trees. The test measurement was carried out on May 23-24, 2018.

ALS Data
In the presented work, data from ALS acquired in 2012 as part of the ISOK was used. which included the area of the study, in the presented work. The obtained data was characterized by the density of points at the level of 12 points per square meter. The data are provided grouped into files, each corresponding to the extent of the sheet of map in the flat rectangular coordinate system "1992" in the scale 1:1250. The area of study was covered by four sheets: The study uses data in the LAS format. The LAS (LASer File Format) format is a public format for handling and exchanging 3D point clouds between users. It was developed in 2003 by ASPRS -the American Society for Photogrammetry and Remote Sensing, as the LAS 1.0 Format Standard. It is the world standard for handling laser scanning data and has been approved by manufacturers of equipment (scanners) and software designed for processing data from LiDAR.

Reference Data
For each point reference model from direct tachometry measurements have been obtained. The measurement was made with the TRIMBLE model M3 total station. The device has the ability to perform non-reflective measurements in Direct Reflex technology, which was a very useful option when aiming at high obstacles located near the zenith, where it was impossible to look through the telescope. However, during a very sunny day, the laser spot was hard to see, therefore it is worth to carry out reflectorless measurements in shady places, during cloudy weather or after dark. In addition, a great convenience was the "infinite" sloths, by which the rotation of the telescope and the instrument, when measuring obstacles on subsequent azimuths was very easy. As a result of tacheometry measurement, a RAW file with horizontal and vertical angles was obtained. On their basis, the azimuth and elevation of each measured point of terrain obstacle was calculated. Based on these values, a polygon representing so called open sky on the sky plot was created.

ALS Data Processing
The main task to do in the presented work was to create terrain obstacle models from ALS data. The input data was four files of point clouds in the LAS format, as mentioned before. For the processing of the ALS data set, a free package created for LiDAR data processing -LAStools was used. This package consists of a dozen or so tools and allows to perform basic operations on point clouds from laser scanning. It was built based on the open source LASlib library (LGPL license) and is available as free for non-commercial purpose. In addition to the LAStools package, another open source software, QGIS, was used to further develop the data. Terrain obstacle model creation was carried out using a dedicated author's QGIS plugin named ALSObstModel.
The obstacle modeling process consisted of several steps. The first of these was to combine all of point clouds covering the area of study into one, using the lasmerge tool of the LAStools package. The resulting cloud was 1 GB in size and was time and sources consuming in further process. Therefore, in the next step, a data set limited by a cylinder was prepared for each considered measuring point. The following assumptions were made during the cutting out of data: -cylinder axis is parallel to the vertical and passes through a measuring point -the radius of the cylinder ensures the inclusion of all the objects higher than 10 degrees above the horizon of the GNSS antenna.
Calculation of the cutting radius took place in two stages. At the beginning it was assumed that there are no objects higher than 50 m in the area under consideration and the initial cut off was made. Next, the height of the highest point in the obtained area was determined and on this basis the definitive radius of cutting (ds) was determined as follows:   The model prepared in this way is saved in the text format and can be used to plan GNSS sessions as well as to pre-process measurement data. Generated files in next step were loaded as CSV text layers into the QGIS application. On their basis, the points were combined into polygons, so that it was possible to compare them with polygons developed from tachometric measurements. The created objects represent an open sky as in tacheometry.

Analysis of Results
From the point of view of determining the GNSS satellites visibility, the measure of the quality of the terrain obstacle model is the accuracy of mapping the boundary between the obstacle and the open sky. In order to evaluate the models obtained from the ALS data, these models were compared with the models obtained from direct tacheometric measurement considered as reference. QGIS software was used to compare the obtained models.
The accuracy analysis was performed by comparing the polygons representing the open sky. The analysis was performed in the QGIS software. First, using the geoprocessing/symmetric difference tool, a polygon layer was created which is the difference of the sum of polygons and their common part. As it can be seen, in case of each of the models there are areas characterized by errors with large values. Comparing the obtained polygons with hemispherical images of the surroundings of the measurement points, these errors can be classified as gross errors. The source of these errors is twofold. The yellow color indicates the area where errors resulted from the fact that the tree was cut, which grew when the ALS data was acquired. In the red rectangle, however, areas where errors resulted from the non-mapping of the building wall from ALS data were marked. The lack of building wall on ALS obstacle model results directly from the geometry of the air scanner measurement. As a result of such location of the scanner in the resulting cloud there is a high density of points on the roofs and relatively small density on the walls [9]. The number of wall points depends on the direction and speed of the flight, as well as the scanning density. However, as it can be noticed in the attached figures, regardless of the direction of the flight, densities of the points on the walls are too low to create a suitable model.
In the next step elevation errors for subsequent azimuths (1 degree resolution) were calculated. Taking into account all data appeared errors reaching 60 degrees, marked on the chart in red. However, if we remove data with gross errors from the data set, we will obtain elevation accuracy at the level of single degrees. The maximum error was 22° at point 0838806, 12° at point 0839104 and 9° at point 0838904, while the average errors were 0.6°, 0.7° and 1.7° respectively.

Conclusions
The paper presents a method of terrain obstacle modelling based on LIDAR data and analysis of the obtained models. Modelling was carried out using open source software. The main advantage of the ALS data is the possibility of modelling for large areas and the lack of a need for a direct visit to the measuring point. The problem is, however, the timeliness of the data. An example of this is the excision of a tree in the vicinity of one of the analyzed points, which caused gross errors in model. The second issue to deal with is difficulty in mapping the walls of buildings with ALS data. The assumption of the proposed method was the use of a raw point cloud and the maximum simplification of the modeling process. In order to cope with the problem of building walls modelling, it would be necessary to add building models from external sources or apply the algorithm of building modelling, which will be carried out in future research.
As a result of the modeling, polar plot of the elevations of the obstacles above horizon have been obtained directly from ALS data. Apart from the abovementioned gross errors, an average error in determining the elevation of the terrain obstacles from 0.6° to 1.7° was obtained. The road equal to one elevation over the horizon is covered by the satellite in about 3 minutes, so the accuracy obtained allows for the planning of satellite visibility with an error of up to 6 minutes, which in the case of planning static GNSS sessions is a permissible quantity.