CREATION OF A MULTIRESOLUTION AND MULTIACCURACY DTM: PROBLEMS AND SOLUTIONS FOR HELI-DEM CASE STUDY

The work is part of “HELI-DEM” (HELvetia-Italy Digital Elevation Model) project, funded by the European Regional Development Fund within the Italy-Switzerland cooperation program. The aim of the project is the creation of a unique DTM for the alpine and subalpine area between Italy (Piedmont, Lombardy) and Switzerland (Ticino and Grisons Cantons); at present, different DTMs, that are in different reference frames and have been obtained with different technologies, accuracies, and resolutions, have been acquired. The final DTM should be correctly georeferenced and produced validating and integrating the data that are available for the project. DTMs are fundamental in hydrogeological studies, especially in alpine areas where hydrogeological risks may exist. Moreover, when an event, like for example a landslide, happens at the border between countries, a unique and integrated DTM which covers the interest area is useful to analyze the scenario. In this sense, HELI-DEM project is helpful. To perform analyses along the borders between countries, transnational geographic information is needed: a transnational DTM can be obtained by merging regional low resolution DTMs. Moreover high resolution local DTMs should be used where they are available. To be merged, low and high resolution DTMs should be in the same three dimensional reference frame, should not present biases and should be consistent in the overlapping areas. Cross-validation between the different DTMs is therefore needed. Two different problems should be solved: the merging of regional, partly overlapping low and medium resolution DTMs into a unique low/medium resolution DTM and the merging with other local high resolution/high accuracy height data. This paper discusses the preliminary processing of the data for the fusion of low and high resolution DTMs in a study-case area within the Lombardy region: Valtellina valley. In this region the Lombardy regional low resolution DTM is available, with a horizontal resolution of 20 meters; in addition a LiDAR DTM with a horizontal resolution of 1 meter, which covers only the main hydrographic basins, is also available. The two DTMs have been transformed into the same reference frame. The cross-validation of the two datasets has been performed comparing the low resolution DTM with the local high resolution DTM. Then, where significant differences are present, GPS survey have been used as external validation. The results are presented. Moreover, a possible strategy for the future fusion of the data, is shortly summarized at the end of the paper.


INTRODUCTION
Elevation data are fundamental in Geographical Information Systems (GIS) and one of the spatial data provided by the European INSPIRE (Infrastructure for Spatial Information in Europe, 2007) Directive.They are stored in Digital Elevation Models (DEM), Li et al., (2005), and El-Sheimy et al., (2005), which are the basic tools for many environmental and territorial applications.DEM refers to the generic family of elevation models: they are distinguished in Digital Surface Models (DSMs), which represent the actual surface (including buildings, woods, etc.) and Digital Terrain Models (DTMs) which represent the elevation of bare soil.DTMs are important in environmental engineering, civil engineering and in many specific applications, such as the design of telecommunication networks or precise agriculture (Pirotti et al., 2013a), especially in the advent of high-resolution sensors for 3D modeling (Pirotti et al., 2013b) which allow accurate discrimination of vegetation from ground plane (Pirotti et al., 2013c).In recent years, both in Italy and Switzerland, elevation data with different reference frames, technologies, accuracies and data models have been acquired.HELI-DEM (HELvetia-Italy Digital Elevation Model) project (Biagi et al., 2011), funded by the European Regional Development Fund within the Italy-Switzerland cooperation program, aims at developing a unique DTM for the alpine and subalpine area between Italy (Piedmont, Lombardy) and Switzerland (Ticino and Grisons Cantons).The DTM will be correctly georeferenced and produced validating and integrating the data that are available for the project.To meet this task with a horizontal resolution equal or better than 20 meters and the best achievable accuracy according to the available products, a recomputation of the geoid, i.e. the elevation datum, with an accuracy better than 10 cm is required (Barzaghi et al., 2007, Marti et al., 2007).Moreover, a permanent GNSS network has been set up to validate and georeference GNSS surveys along the borders.The institutions involved in the project are: Fondazione Politecnico di Milano, Politecnico di Milano, Politecnico di Torino, Provincia Autonoma di Bolzano, Regione Lombardia, Regione Piemonte and SUPSI (Scuola Universitaria della Svizzera Italiana).Two other institutions, SwissTopo and the Istituto Geografico Militare Italiano (IGM) are involved as external partners.

Available data
The first aim of HELI-DEM was the collection of data and metadata for the study area.In particular, this phase involved the research and the collection of local DTMs, i.e. digital terrain models with medium to high resolutions, developed by local authorities (Regions, National Geodetic Authorities), with a national or regional extension.The following data have been made available to the project partners:  official DTMs of Lombardy, Piedmont and Switzerland with resolutions of 20, 5 and 25 meters respectively,  DTM with a 10 meters horizontal resolution for the eastern area of Piedmont,  Lidar DTM (Ackermann F., 1996) that covers the catchement area of the Como Lake with a horizontal resolution equal to 2 meters,  high resolution Lidar DTM (1 meter), produced by the Ministry of Environment, which covers the Lombardy and Piedmont main hydrographic basins.The DTMs are in different reference frames.The Lombardy DTM is in the Italian Roma 40 reference frame, projected in the Gauss Boaga (West zone) cartographic coordinates (R40-GB), while the other DTMs are provided in ETRF89-IGM95.In particular, the two DTMs of Piedmont and the Lidar DTM of Como Lake are in UTM coordinates, while the others (the Swiss DTM and the HR Lidar DTM) are in geographic coordinates.Figure 1 shows the distribution of the collected DTMs over the project area of interest.The creation of a unified multiresolution DTM for the alpine and subalpine area between Italy and Switzerland is the main goal of HELI-DEM.Therefore some conclusions about the final resolution and accuracy of the resulting model can be drawn depending on the available DTMs.To this purpose, (Biagi et al., 2012a) shows a graphical representation of the maximum resolution available for the provinces involved in the project.In order to create a correct multiresolution DTM, several cross-validations have been carried out.

DTM
To cross validate the DTMs, two different analyses have been done.The first is the comparison between DTMs with similar horizontal resolution where they overlap, i.e. in the borders between their domain.The second comparison regards low and high resolution DTMs where the second one exist.The results of the comparison between low resolution DTMs in the overlap areas have been already presented in (Biagi ed al., 2012a).
Here we just mention that they do not present biases and the general statistics are consistent with the nominal accuracies of the datasets.The purpose of this paper is the cross-validation of low resolutions (LR) and high resolution (HR) DTMs Lombardia region with a particular detail to the case study Valtellina valley.Moreover, were differences are anomalous an external cross-check will be performed by RTK GPS (###).At the present, we are merging low resolution DTMs to create a unique DTM for the alpine and subalpine area between Italy and Switzerland.The results will be presented in a following paper.

VALIDATION OF HIGH RESOLUTION DTMS IN LOMBARDY
HR and LR DTMs have been compared by implementing new commands in GRASS 7.0.The procedures that have been applied for the cross-validation in the Lombardy region are here shortly described.
The LR DTM has been transformed to the reference frame of the HR DTM (ETRF89).The result is a list of 3D points, that clearly are almost regularly spaced but no more on a regular grid.For each point of this dataset a square window of 20 × 20 m 2 has been considered.All the nodes of the HR DTM contained in the window are averaged and their mean is compared with the LR point.Note that the HR DTM covers only the hydrographic main basins of the region, for a total coverage of about 3810 Km 2 (which corresponds to 16% of the total Lombardy region).
The following statistics have been obtained: = -138 m To verify the existence of global biases, the differences have been clustered in six classes: the results are reported in Table 2.No global bias is present and general statistics are consistent with nominal accuracies.Outliers are present but in negligible percentages.Figure 1 depicts the spatial distribution of the differences.Moreover we have analyzed a phenomenon that, at first sight, seems quite strange.As example, Figure 2 shows the differences in the Valley of S. Giacomo that covers all the B1 Lombardy grid.The mean difference is 0.5 meters, that is not significant with respect to the LR DTM accuracy.However, in the West slope (the lower part of the figure) a negative bias between HR and LR DTMs appears.In the other slope a positive bias is present.Such a behavior could be justified by a translation of the two DTMs with respect to the valley axis.A program to estimate biases between DTMs has been implemented to verify the existence of such translation.

Class
Limits Percentage  3. HR-LR height differences in the valley of S. Giacomo  The estimate of a translation and a bias between two DTMs is performed by Least Squares (Koch, 1987).The height of a point P is described by the following: where 12 , ff are respectively the height functions of DTM1 (in our case, the HR DTM) and DTM2 (the LR DTM), υ 1 , υ It can be linearized by the following where f  is the gradient of the function.By all the observed differences , the translation and the bias can be estimated by LS.To fill the design matrix, the gradient of the LR DTM are computed by the usual numerical approximation.The translation and the height bias have been estimated in two different ways: firstly, considering the Lombardy official cartographic partition in sheets and then considering the natural orographic partition in valleys.In the first case the estimates for the individual sheets are always smaller than ±1.5 m, and their mean is 0.3 m.No systematic translation or bias appear at the regional scale or in the individual sheets.Also in the case of the single valleys, no significant translations or biases exist, with estimates alwas smaller than 4 and 2 meters respectively.These results are acceptable if compared with the nominal accuracy of the regional LR DTM.Therefore, no translations or biases exist, neither at the region scale nor for particular valleys.However, local biases seem to be present, where the differences of LR with respect to HR DTMs are significant and spatially correlated.Therefore, a deeper check has been done, by two different approaches: 1. at first, DTMs have been compared with the monographs of the zero order Italian geodetic network, 2. then, RTK GNSS surveys have been performed in areas where HR and LR DTMs differences exceed 10 meters.
The tests have been performed on a case study, that is Valtellina valley (Figure 1).The choice of this area as study-case is due to specific characteristics: Valtellina is the most important and most populated Lombardia valley in the alpine area and, furthermore, it is localized in the northern area of Lombardia region, close to the Swiss border, the perfect area for our studies.

DTM comparison with IGM95
IGM95 is the zero order geodetic network, designed and surveyed from the '90s by Istituto Militare Geografico Italiano (IGM, the official cartographic Italian authority, Surace, 1997).Now, the network consists of about 2000 benchmarks, that are typically monumented on the ground or on artifacts.We got a set of IGM95 benchmarks in the study area.We have carried out two comparisons on their official monographs: the former with the low resolution DTM, the latter with the high resolution (LiDAR) DTM.Only recently re-surveyed IGM95 benchmarks have been used for the check, in order to guarantee their cm level accuracy: in total 64 points are available, of which 6 in Valtellina.In the comparison of LR DTM, the mean difference is 1 meter and the standard deviation is 3 meters: these results are consistent with the accuracy of the DTM.The comparison with HR DTM requires more detailed analyses.The mean difference is of 6 cm, with a standard deviation of 13 cm.Just 13 points present differences greater than 30 cm: all of them are monumented on artifacts that clearly are not in the DTM, like, for example, walls or road embankments.On this regard, despite the small number of check points, these results are completely satisfactory.However, no IGM95 benchmarks are present in the area where big differences exist between LR and HR DTMs; therefore, ad hoc RTK GNSS surveys have been planned and that will be described in the following section.

DTM comparison with RTK surveys
We want to confirm the accuracy of HR DTM: for this reason GPS surveys can be used as external validation in selected areas.RTK (Real Time Kinematic) surveys seem adequate, because they provide accuracies typically better than one decimeter.A permanent GNSS network of reference stations has been created ad hoc for HELI-DEM project, by selecting a subset of already existing stations in Italy and Switzerland.After careful analysis of the state of art, almost all the permanent stations located in the area of the project and some stations present in neighboring areas (for reasons of geometric configuration) have been chosen; arrangements have been made with the owners of the stations (GPSLombardia, GNSSPiemonte, GeoTop and AGNES) to obtain the real time data flow: more details on the HELI-DEM permanent network, its purposes and characteristics can be found in (Biagi et al., 2012b).After comparing LR and HR DTMs, 12 study areas have been chosen, where height differences are bigger than 10 meters and are significantly correlated in space (Figure 1).The area identification has been performed also by using Google Earth in their pre selection.In particular, the study areas have been chosen also on the basis of a good access and potential GNSS visibility criteria.Moreover, we have chosen areas where the ground seemed smooth, to avoid errors due to a bad definition of the medium ground.Areas VT01-VT02-VT03-VT04 are localized in west side of the valley, between Colico and Sorico, while areas from VT05 to VT12 are in east side between Sondrio and Bormio (Table 4).VT01 to VT06 are below 500 m of elevation, VT07 to VT10 are from 500 m to 1000 m., VT11 and VT12 are above 1000 m.RTK surveys have been made with Leica RX 1200 receivers.Rover receivers communicated with HELI-DEM permanent network by UMTS connections.Contrary to the surveying practice, we did not monument benchmarks but just surveyed height profiles on given paths.Indeed, for each study area, we are interested to check heights on the whole area and not on individual points.We moved along roads or paths easily accessible.Along flat sections, fast static measures of 5 seconds were repeated every about 20 meters.In slope sections, measures were more frequent, at about 10 meters.Each path has been surveyed twice, in order to have a cross check of the results.Moreover, all the surveys have been repeated in two campaigns, the former in June/July, the latter in October 2012.Survey conditions were not always optimal, due to the presence of vegetation and obstructions.Therefore, in some places, ambiguity fixing was not possible.In any case, for each study area, at least 50 points have been surveyed in each campaign: in total, 1300 RTK points have been collected.Considering the alpine nature of the case study, these results are completely satisfactory.Colico zone (areas V01, V02 and V03) provides the only one exception.Indeed, in October campaign, few points were collected due to significant problems in GNSS data acquisition in the whole zone.Probably, this was caused by electromagnetic disturbance emitted by some repeater.Data collected during surveys, stored in memory cards, were downloaded and analyzed using Leica Geo Office software.The software allows to visualize, import, manage and process GNSS data and RTK results.Then, RTK results were converted to text files for further analyses with other GIS and statistic software.It is worth to note that RTK results are in ITRF2008 while HR DTM is in ETRF89.Before comparing them, RTK results have been transformed to ETRF89 by the standard transformations (Donatelli et al., 2002, Boucher andAltamimi, 2011).

Id
After the transformation, HR DTM and RTK points have been compared in the following way.Each RTK result consists of the 3D coordinates of a point measured on the ground.A simple program has been implemented: its inputs are the DTM grids and a list of the RTK estimated coordinates.For each RTK point, the program searches the nearest sixteen nodes of the DTM.The heights of the nodes are interpolated with a bicubic exact interpolation in the RTK horizontal coordinates.RTK height is transformed from geodetic to orthometric by applying the ondulation of the official Italian geoid Italgeo08 (Barzaghi et al., 2007) Therefore, the difference between interpolated DTM and RTK is computed.Clearly, in this way, horizontal RTK errors are neglected: however, for the purposes of this comparison, they can be considered as not significant.About ten differences bigger than 1 meter exist.In any case we consider all the 43 height differences greater than 50 cm as outliers and pay particular attention to them.  5 and Table 6. Figure 6 provides the example of VT09.

MERGING OF LR AND HR DTMs
To achieve the goal of the project, two different DTMs have to be produced: both of them will have a horizontal resolution of about 20 m (therefore, can be defined Low Resolution) and will cover the whole area of the project.Note that, in order to optimize the numerical computation implied by steps 4 and 5, a Fast Fourier Transform (Oran Brigham, 1974).approach will be applied.Firstly, FFT of the differences will be computed, then they will be multiplied by the Fourier transform of the filter, at the end the inverse FFT of the product will be computed.In that way we will obtain the LR DTM corrected using the HR information coming from the LiDAR DTM.

CONCLUSIONS
The presented research is in the framework of the Interreg Italy-Switzerland HELI-DEM project, that is aimed at creating a transalpine DTM by merging all the available height information for the alpine area.The paper discussed the first operations needed for the final merging of HR and LR DTMs.Several cross-checks have been performed: firstly, LR DTMs have been compared in their overlapping boundaries.Then, they have been compared with HR DTMs where they are available.In the paper, the cross check of LR and HR data for a case study in Lombardy region (Valtellina valley) has been presented.In some areas, Lombardy LR and HR DTMs show an apparent horizontal translation.A software to estimate the translation between two DTMs by least squares has been implemented.Globally, no common translation exists: this is also true at the level of the several subareas in which the valley can be split.In order to have an external check, the HR DTM has been compared with the published coordinates of the zero order Italian geodetic network: the results are completely satisfactory.More analyses have been carried out in area where the differences between the LR DTM and the HR DTM are bigger than 10 meters, by performing several RTK GNSS surveys.Except for two areas, that are discussed in the paper, these gave completely satisfactory results, that are consistent with the nominal accuracy of the HR LiDAR DTM.The future work will consist in the merging of all the input data into a common final DTM: this will be reported in following papers.

Figure 1 .
Figure 1.Map of the height differences between HR and LR DTMs in San Giacomo Valley.

Figure 2 .
Figure 2. HR-LR differences in the valley of S.Giacomo.The presence of systematic patterns is clear (biases with opposite signs on the valley left and right slopes are hardly observable).
A first, visual, inspection of RTK/DTM profiles evidences that most outliers are clearly caused by RTK survey problems and do not represent DTM blunders.As expected, typically they are present under dense vegetation or in obstructed sites.Part of them regards sparse points, the others affect close points in small, particularly obstructed, areas and repeat in both the campaigns.Two exception exist.VT07 (Errore.L'origine riferimento non è stata trovata.) is in Lovero (SO): the path is about 900 m long, with a total height difference of 45 m.During the first campaign, 7 close RTK results (brown and red in figure) are clearly blundered while all the other Δh are between 0 and 30 cm.This problem does not repeat in the second campaign.It is clearly due to a prolonged false ambiguity fixing during the first campaign that happened in an area that, from a GNSS observational point of view, should be defined good.VT10 (Figure5) is in Sondalo (SO): a paved road 600m long that climbs orthogonally a steep slope of the valley (total height difference of 44 m).RTK heights are systematically higher than DTM heights in the road side toward the valley and this survey contains most of the differences bigger than 1 m: this bias is clearly caused by a small error in the DTM horizontal georeferencing.The results of the other differences are satisfactory and are shown in Table

Figure 4 .
Figure 4. VT07 area: GNSS RTK results and differences with respect to HR DTM heights.

Figure 5 .
Figure 5. VT10 area: GNSS RTK results and differences with respect to HR DTM heights.

Table 1 .
The available data.A: nominal vertical accuracy.S: spatial resolution.C: covered area within the HELI-DEM region

Table 4 .
Positions of the RTK study-areas.L: length of the path, H: orthometric height, : mean ITRF2008 latitude and longitude of the area,

Table 5 .
Statistics of first campaign.Id, Number of RTK Points, Number of Outliers, Percentage of Outliers, Mean, Standard deviation, Absolute Maximum value, Root Mean Square Error.All metric results in cm.

Table 6 .
Statistics of first campaign.Id, Number of RTK Points, Number of Outliers, Percentage of Outliers, Mean, Standard deviation, Absolute Maximum value, Root Mean Square Error.All metric results in cm.