The LiDAR-based topo-hydrological modelling of the Nigula mire, SW Estonia

There was no big influence of the used cell scales and algorithms on the mean topographical properties of the Nigula mire digital elevation models (DEMs). The DEMs, generated using the Triangulated Irregular Network and Inverse Distance Weighted algorithms, revealed the closest mire surface properties from all used generation algorithms. The subtracted MAXMIN DEMs layer revealed a well visible net of ditches and possible plant cover pattern differentiated in vertical scale. In the Nigula mire 58% of the mire surface basins have SSW orientation, followed by the levelled and less fractionated NNE basin region (23% coverage), and the most fractionated but with steeper sloping W-orientated basin region (8% coverage).


INTRODUCTION
Ground and surface water, in its specific dynamics, has a major role in the development of mire landscapes and their surface patterning. Already in the early 20th century it was reported that the vegetational patterning of raised bog was closely related to landform morphology and hydrology, advancing so the development of a threedimensional and dynamic view of the raised bog development (Glaser 2002). Since then, the general term of patterned peatland or patterned mire has been used for sites with a specific type of pattern, formed by mire microforms for instance. They can be found on bogs, fens or a combination of the two (Rydin & Jeglum 2006).
Peatlands patterning is closely related to surface water discharge, not only across the mire boundary into the surrounding landscape, but also to the processes occurring within the mire which regulate the flow towards the boundary (Ingram 1983). In the case of elongated features of the mire surface pattern, usually stretched out perpendicular to the mire surface slopes, these processes can provide useful empirical information about the direction of the main flow paths in the mire area (e.g. Ivanov 1981;Rydin & Jeglum 2006). Price & Maloney (1994) concluded that the movement of water within, and out of, patterned peatlands is strongly controlled by the nature and position of pools and ridges within the peatland water basin. In a patterned bog with relatively high ridges the ridges can act as local aquitards in the mire, hydrologically separating pools from each other and also causing the expansion of the depressions and storage of surface water retention in the area.
Changes in mire microtopography and in situ distribution of mire surface microforms or complexes of those (i.e. microtopes), together with the corresponding surface roughness properties, have been shown to be important hydro-topological factors in the formation of discharge components for the peatland catchments. The surface runoff or overland flow, seepage, groundwater and pipe flow, and open channel flow (i.e. rills, streams or even rivers), have been named as the key hydrological pathways of those catchments (e.g. Eggelsmann et al. 1993;Holden 2005).
However, because of patterned and at the same time comparably levelled surfaces the delineation of dischargeforming catchments on the peatlands is difficult (Ingram 1983;Eggelsmann et al. 1993). Adding for all, a mire macro-landscape (Galkina 1946, cited in Masing 1998 or a system of mire massifs (Ivanov 1981), or mire complexes (Moen 2002) surrounded by minerotrophic fens of variable extent, could consist of several elementary mire massifs (Galkina 1946, cited in Masing 1998, where the residual and secondary developed water bodies, lagg fen areas and mineral islands exist. Hydro-logically these mire massifs in complexes are connected via different mire types and may contain Carex fens, wet forests, flat bogs and string mixed mires (Rydin & Jeglum 2006). Very often the general shape of a mire complex and its surface patterning also reflects interactions between the underlying mire terrain form (e.g. paludified mires), and the regional climate and hydrological environment of the peatland location (e.g. Ingram 1983).
About two thirds of Estonian mires have developed due to paludification of mineral land (incl. 30% of mire formation on sandy soils) (Allikvee & Ilomets 1995), and about one third have formed by 'infilling' of the surface water bodies, also called terrestrialization or hydroseral succession (e.g. Lode 1999). According to the general surface shape classification of elementary mire massifs, in Estonia there are mires with (1) flat or slightly convex sloping surfaces, sometimes located on the sites of spring upwellings, (2) flat or gently undulating surfaces of fens or transitional mires, i.e. poor fens, and (3) convex surfaces of bogs, which can be divided into (a) slightly convex surfaces with minimally segregated marginal slopes, (b) flat plateau surfaces with deep marginal slopes and (c) strongly convex surfaces with extensive and obvious marginal slopes (Masing 1988).
Currently, delineation of wetland basins and evaluation of water discharge dynamics are an important tool for integration of wetlands (incl. peatlands or mires) into the EU Water Framework Directive (EC 2003).
The LiDAR (short for Laser Imaging Detection And Ranging) data sets existing in Estonia provide highquality three-dimensional surfaces of patterned peatlands and their delineated water basins. The availability of elevation data in digital format has favoured development of automated tools that can be used to delineate drainage basins and their associated stream network (Furnans & Olivera 2000). Accurate digital surface data, visualized via digital elevation models (DEMs), are the primary advance of LiDAR data, which are less subject to the horizontal error inherent in comparison of contour lines derived from ordinary contour maps (Haile 2005). Krause & Bronstert (2005) expressed their doubts about the possibility of using algorithms of automatic digital terrain analysis for delineation of watersheds within lowland floodplains, because of a minimal extent of topographical heterogeneity required for the generation of DEMs. However, analysis of laser scanning of the fine-scale pattern along hydrological gradients in a peatland ecosystem has given promising results in a range of ecosystems where the vegetation pattern is linked to the ecological function (Anderson et al. 2010).
Many factors affect the accuracy of DEMs, i.e. the accuracy, density and distribution of the source data, the algorithms used for data interpolation and the DEM resolution (Liu et al. 2007;Liu 2008). Anderson et al. (2006) indicated that DEM horizontal resolution significantly influenced the level of reduction that LiDAR data sets could withstand, although for the satisfied soil-landscape models the DEM resolution could vary from 2 m resolution (Gessler et al. 2000) to 10 m and 30 m resolutions (Thompson et al. 2001). Results of the LiDAR DEM modelled flow network in Murphy et al. (2008) showed that the most accurate representation of the actual field-mapped network was even more accurate than the aerial photo-interpreted hydrographical data.
On the other hand, results of LiDAR-derived DEM modelling of the most important hydrological features, such as drainage network and boundaries of the subbasin (also called basin or watershed), showed high sensitivity to both DEMs accuracy and the resolution used (Liu et al. 2005). However, it was demonstrated that the LiDAR-derived DEMs of high accuracy and high resolution offered the possibility of improving the quality of hydrological features derived from DEMs (Ibid).
The aims of this paper are (a) to compare the topographical properties of different LiDAR-derived DEMs (generated for different pixel sizes with different algorithms) for the Estonian Nigula mire complex, (b) to demonstrate the mire basin modelling results of different LiDAR-derived DEMs, (c) to compare the LiDARderived DEM modelling results with results derived from an empirical map at a scale of 1 : 10 000, constructed from results of earlier field measurement campaigns (presented in Raukas & Kink 1993/94) and (d) to define the most suitable topo-hydrological modelling approach for the Nigula mire landscape analysed.

Study area
The Nigula mire is located on the Pärnu lowland, in the southwest of Estonia, 10 km to the east of the Gulf of Riga (Fig. 1). The peculiarity of the surface pattering of the contemporary Nigula mire is related to the origin of the mire formation. It began as an 'infilling' of Ancylus Lake during the Boreal period, first in the western part and thereafter in the eastern part of the lake (Pirrus 1963;Karmu 1966;Loopmann et al. 1988). About 8000-8500 years BP the western part of Ancylus Lake was 'infilled' by the mire. The terrestrialization of the eastern part of the lake started about 5000 years later, i.e. 3000-3500 years BP (Raukas & Kink 1993/94). The NE-SW orientation of the current Nigula mire is defined by the orientation of the small-drumlin thresholds from the Ancylus Lake period, at present partly 'buried' under the mire, partly outcropping as small-drumlin residuals, forming mire mineral isles. Five outcropping mire isles occur in the current Nigula mire area, of which four, having the sequential NE-SW orientation, divide the mire into two parts -one third of the area is the western Urissaare bog and two thirds of the area make up the eastern Nigula bog (Kolla 1982) (Fig. 1).
A surface track depression of the previous Lemme River bed, beginning on the pool-ridge microtope in the middle of the Nigula bog, is located between the mineral island to the east and Lake Nigula to the west. The depression runs towards the SW corner of the Nigula mire. According to a field survey carried out in the middle of the 20th century, the Lemme River bed depression was hardly segregated on the bog surface (Karmu 1966), but in the undercut stream bed region the breadth of the stream bed was about 1.0 m and the depth about 0.5 m (Puura et al. 1990). The velocity of the open water river flow at its start was about 0.8 m/s, but increased towards the SW corner of the mire area up to 2.3-3.0 m/s (Karmu 1966).
The Nigula bog is a typical Estonian southwestern large bog-type bog with a flat plateau. The bog has a relatively steep 0.03-0.05‰ marginal slope in the east; the other marginal slopes vary slightly around 0.002‰ (Loopmann et al. 1988). The average slope of the southern oriented plateau is about 0.001-0.0015‰ (Karmu 1966). The western marginal slope of the Urissaare bog, contrary to that of the Nigula mire, is comparatively higher (2.3 m) and steeper (0.03-0.05‰), but both these values decrease towards the south. In the eastern part the bog is smoothly connected to the mineral island (Karmu 1966;Loopmann 1970) (Fig. 1).
The highest point of the Nigula mire is located in the northern part of the Urissaare bog (59.5 m a.s.l.) (Loopmann 1970). The average peat depth over the whole Nigula mire is about 3-5 m with a maximum of 6-7 m (Karmu 1966).
The prevailing mire type in the Nigula mire is the bog. Only mineral isles and narrow belts of the transitional mire skirt the dominating bog-type mire areas (Puura et al. 1990). However, two non-bog type parts exist in the Nigula bog area: (1) the NE 'corner' of the mire, called a quaking or hag-type mire (Loopmann 1970), and (2) the SE 'corner', where probably some hand peat-cutting took place at the beginning of the 20th century (Karmu 1966).
The main area of the present mire surface is fed by precipitation. The excess water from the overland flow is collected in the belt areas between the mineral isles and the bog massifs, and upwelling of the bog karstik springs occurs not far to the east of the northern area of the mineral isles (Karmu 1966).
The Nigula mire is a water divide into four main surface water courses, called the Salatsi, Lemme, Häädemeeste and Rannametsa rivers. The Lemme and Häädemeeste rivers run out of the Nigula mire towards the SW, the Rannametsa River flows towards the NE and the Salatsi River in a SE direction (Karmu 1966). According to Loopmann (1970), the total length of the ditches surrounding the Nigula mire is about 17-25 km, the average breadth of which is 2 m and depth 1 m. The fluvial surface water runs into the ditches mainly during the spring and summer-autumn flooding periods.
The GIS-based study area, considered in this paper, forms 3845.10 ha (perimeter 26.08 km), covering the 2411.03 ha Nigula mire area (perimeter 35.34 km). The GIS study area is limited by roads to the north, east and south, whereas the eastern road follows the course of the Tuuliku drumlin ridge. In this study the Nigula mire itself, as a landscape, is defined by the '0' contour line of the peat soils (Fig. 1).
The Nigula mire pools, with open water surfaces, have been classified as hillside and water-divide pools, pool-lakes, depression and upwelling well pools, and river bed remnant and after-fire pools (Karmu 1966). Currently the Nigula mire is covered with 1261 open water bodies of different sizes (i.e. Lake Nigula, pools and shallow inundated hollows) which were digitized manually from the orthophoto of the year 2005 at the beginning of the current study. The open water bodies cover a total of 71.78 ha (i.e. 3.0%) of the mire area, including the Lake Nigula area of 20.13 ha.
In the GIS-based study area the total length of fluvial water pathways, including ditches, is about 120 km, of which about 50 km (i.e. 41%) is located inside the mire border area, i.e. directly on the mire. The elevation difference between the highest Northwest section (60.83 m a.s.l.) and the lowest Southwest section (49.35 m a.s.l.) of the mire is 11.5 m. The longest longitudinal profile is about 9.5 km and the average breadth in the central part of the mire area is about 3.5 km (Fig. 1).

Creation of DEMs
The LiDAR data were used for different generated DEMs. The laser altimetry data used were gained from a Leica ALS50-II scanner device during a flight survey on 10 May 2008. The average flying height during the survey was 1 km and the resulting First Echo point density at nadir was expected to be 2.3 points/m 2 . After the data post-processing the vertical accuracy of the laser data points was determined to be 8 cm. The data point 'cloud' was then split into 500 × 500 m subsets and classified into Overlap and Erratic High and Low point classes with the TerraScan software. LastEcho points were used to generate the Ground point class, which was the basis for the generation of DEMs (Fig. 2). The ground point data density was determined to be 0.8 points/m 2 . Using the Triangulated Irregular Network (TIN) and Inverse Distance Weighted (IDW) algorithm implementation in the LP360 module (i.e. the LiDAR extension for ArcGIS, RockWare Inc.), several DEM surfaces with different pixel sizes (1 × 1 m, 2 × 2 m, 5 × 5 m, 10 × 10 m) were generated (later named as TIN or IDW DEMs).
The Local Binning Algorithm (Kim et al. 2006) with a standard Search radius (equals cell resolution 2 2 m, × in our case the suggested cell resolution is 0.63 m), was used for the creation of MAX and MIN DEMs, whereby the maximum values track the highest elevation points, usually to the top of vegetation, and minimum values to the ground level surface (Ibid). Since all altimetry data were collected in the GRS-80 height system, the surfaces were then normalized to a Bk77 height reference system, using a geoid undulation of 19.75 m.
A total of ten LiDAR-based DEMs were generated for the following 3D surface analysis and Hydrological Basin modelling.
Empirical (EMP) DEMs were generated with the ArcMap software via the TopoToRaster module, whereas the elevation points for modelling were derived from a digitized empirical map layer at a scale of 1 : 10 000, constructed from results of earlier field campaigns and previously presented in Raukas & Kink (1993/94). In total four EMP DEMs of 5 m and 10 m cell resolution were created with differently defined Ditch Feature Layers, e.g. with and without SW ditches (Fig. 1). The used density of the Contours was 0.2 m in generation of EMP DEMs. The model network layers of Streams and Roads were extracted from the Estonian Base map and were elaborated with corresponding visualizations from the orthophoto of 2005 (Fig. 1). The Lake layer was obtained from manually digitized water bodies of the orthophoto 2005 layer (Lode et al. 2011).
All DEMs were created within the GIS-based Study area covering the whole Nigula mire landscape. A buffering belt was left between the net of the main road system and the mire '0' contour line (Fig. 1). MIN and MAX DEMs with 10 m cell resolution only were used both for the Sink Depth analyses and the Basin delineation runs.

Topographical characteristics of DEMs
All generated DEMs were visualized via the HillShade module. The HillShade layers of visualized 1 m, 2 m, 5 m and 10 m cell resolution DEMs revealed the cell size resolution-dependent mire surfaces, where the 1 m and 2 m cell resolution layers fitted well for visualization of mire micro-topography, and the 5 m and 10 m cell resolution layers for the whole mire landscape (Fig. 3), where manually digitized open water bodies from the orthophoto accord well to the microtopographical structure of the visualized HillShade layers.
There was no noticeable cell scale dependence or generation algorithm impact on General Statistics of the modelled Elevation and Aspect layers for different DEMs over the whole mire landscape (Table 1) Similarly, no differences were observed in the visualized Elevation layers of the different DEMs, contrary to the modelled Slope layers. Almost none of the LiDAR-based Aspect layers had well-observed aspect orientations in the visualized layers. Only the Reclassified Aspect layer to the five classes of 10 m cell resolution DEMs resulted in a visibly meaningful mire surface orientation pattern.
On average, the EMP DEM Elevation and Slope layers showed the lowest surface property values in comparison with LiDAR-derived DEMs (Tables 1 and 2). The surface mean Aspect of the EMP DEM was the most eastern oriented and the visualized Aspect layer had the coarsest, but most easily observed surface cellaspect orientations of all the modelled DEMs (Fig. 4).
Although there was no clear differentiation in general statistics of the EMP DEM from altimetrygenerated DEMs of 10 m cell resolution (Fig. 5), differences were observed in the tree cluster results (Fig. 6), where the Linkage Distance d was shortest between the TIN and IDW DEM data sets (d = 19), and the longest between data sets of LiDAR-derived MAX DEM and the EMP DEM (d = 209 in Table 3).

Basins of delineated DEMs
Results of the Basin delineation of Nigula mire EMP DEMs showed a significant dependence on the ditch features used in modelling runs (Fig. 7). Due to different inputs used for the Ditch Feature Layer in two basin modelling runs, the areal difference for the southern Nigula mire basins constitutes 5.6%, 44.5% and 62.8% (Table 4 and Fig. 7). At the same time, the sum of areal coverage of those basins had differed only by 3.5%, i.e. 62.08 ha. Study of Sink depths of TIN DEMs at cell resolutions of 5 m and 10 m showed that the main areas with high values (e.g. Sink Depths equals 2.04 m for the TIN DEMs at cell resolution of 5 m, and equals 1.38 m at 10 m cell resolution) were concentrated on the main ditch locations. The mean values were found at the sites with tree or shrub coverage and the lowest values were found in the open mire areas (Fig. 8).
Basin model runs of LiDAR-derived 10 m cell resolution TIN DEMs revealed increment of the delineated basin count from 34 to 95 by the change Table 4. Results of two EMP DEM Basin modelling runs for the southern part of the Nigula mire. Depending on Ditch Feature Layer inputs used, in the Basin 1 layer SW ditches were included and in the Basin 2 layer not included in the model run (see also Fig. 7). The cell resolution of input EMP DEMs was 5 m, and 10 m cell resolution was used in the output of delineated Basin layers of Z(limit) from Z(limit) equals max to Z(limit) equals 25% of the max value, and for the 5 m cell resolution from 31 to 55 by the same Z(limit) value changes (Fig. 9). On the other hand, the counts of the modelled basins were rather similar for the TIN DEMs at cell resolutions of 5 m and 10 m (31 and 34 respectively) by the used Z(limit) equals the max value. In the case where the Z(limit) equals 25% of the max value the differences between modelled basin counts at cell resolutions of 5 m and 10 m reached 42% (i.e. 40 basins). Although differently modelled basin layers had different basin counts and basin shapes, there was one main (the largest) basin, delineated from the 10 m cell resolution TIN DEM Basin run, which could be easy to aggregate to almost the same shape and extension, similarly to the EMP DEM Basin run results (Fig. 7), by summarizing the catchments from the other Basin runs (Fig. 9).
The delineated basins for the MIN DEM 10 m cell resolution Basin run were visually similar to the basin delineation at a 5 m cell resolution TIN DEM, and for the MAX DEM 10 m cell resolution Basin run, to the 10 m cell resolution TIN DEM Basin run (Figs 8  and 10). This indicates the better 'reflection' of the finer-scale surface patterning in the MIN DEM layer (although with a lower 10 m cell resolution) in comparison with the MAX DEM layer. All these described Basin runs were performed with max or Default Z(limit) in the Surface Filling procedure before the Basin runs.
It could be stated that there were four main basin regions which could be defined from all Basin runs by the Z(limit) equals Default values for all used DEMs: (1) basin region with a north to northeastern orientation, with a total drainage area of about 560 ha (the minimum number of modelled sub-catchments is 3); (2) basin region with a south to southwestern orientation, with a total drainage area of about 1400 ha (the minimum number of modelled sub-catchments is 11); (3) basin region with an east to southeastern orientation, with a total drainage area of about 270 ha (minimum number of modelled sub-catchments is 4); and (4) basin region with a western orientation, with a total drainage area about 180 ha (minimum number of modelled subcatchments is 17) (Fig. 11). From all the main basin regions the basins with S-SW orientation cover about 58% of the total mire surface area, followed by the N-NE region basins with 23% coverage, E-SE region basins with 11% coverage and W-orientated basins with 8% coverage. The most fractionated basin region is the W region and the least fractionated is the N-NE region.

CONCLUSIONS
The LiDAR data application to the patterned mire landscape in Estonia demonstrated both advantages and challenges of the study of different DEMs, generated directly from the LiDAR raw data material. Not all modelling output results, concerning the influence of different DEMs that used cell sizes and algorithms on different mire topographical and hydrological properties, were presented in the current study, and no ground-level validation was carried out. However, there is a possibility of outlining some general trends covering the research objectives formulated at the beginning of the study.
It has been shown that there were no significant differences in mean values calculated for Elevation and Aspect layers of different LiDAR-derived DEMs for the whole Nigula mire landscape. Thus, in the landscape scale, any of the DEMs could be used for general topographical characterization of the mire surface. However, from a practical point of view (e.g. time and costs), the lower 5 m or 10 m resolution DEMs would be preferred in the case of the Nigula mire.
The HillShade layers, used for hypothetical illumination (ESRI 2007) of generated DEMs, demonstrated clearly the cell resolution-dependent mire surface results, where the DEMs of 1 m and 2 m cell resolution have a great perspective for fulfilling the traditional microstructural study demands on the mire area (Lode et al. 2011). Afterwards it could be concluded that the generated HillShade layers could also be taken as an indicator of the possible cell size scale limits in water basin delineation for the entire mire landscape.
Cell-based evaluation of all generated DEMs depicted higher similarities for the LiDAR-based TIN and IDW DEMs, whereas the generated EMP DEM had the smallest Euclidean distance with MIN DEM. This indicates that MIN DEM can be considered the most suitable DEM generation from LiDAR data sets for substitution of 'ground level' EMP DEM. At the same time the largest Euclidean distance of EMP DEMs from all LiDAR-based DEMs indicates that the empirical elevation map of the Nigula mire reflected less the mire micro-topographical pattern, since the zero level of the mire surface was equalled to the heights of the hummock microform foots (Nastavlenie… 1972) during the field levelling campaign (Raukas & Kink 1993/94).
On the basis of all TIN DEM Basin runs it could be concluded that, at the landscape-scale level, the delineation and count of the modelled mire surface basins are cell resolution and pour point filling property (i.e. Z(limit)) sensitive. Although the Filling of Sinks for the Basin runs is used for removing a '…small imperfection...' in the surface raster data (in the Help section of ArcMap), in the LiDAR data case used, the differentiated Sink Depths, i.e. Z(limit) equals 25%, 50% or 75% of Sink Depth (max), could be regarded instead as modelling options for exercise scenarios of different mire hydrological conditions. In this case, for example, the Basin layers generated by the Z(limit) equals Default (i.e. equals the maximum value of the Sink Depths of the DEM), can be seen as a basin delineation scenario at the highest groundwater level conditions on the mire landscape. But such an approach relies on the assumption of correct raster DEMs generated from LiDAR data.
Subtracted MAX and MIN DEMs of 10 m cell resolution resulted in a well distinguished net of ditches and areal scattered pattering with a height amplitude between 0.00 and 3.54 m, indicating a possible relationship to tree and other higher plant cover distribution in the area with high values and to open areas with the plant cover from the field layer with low values. These results seem to be promising for GISbased identification of plant cover distribution over the mire landscape in the vertical scale.
Two Basin modelling scenario runs displayd delineated mire basins with water divide location on the western border line of the Nigula bog and between the Lemme River bed depression and the same mire western border line. The basins were revealed from the Basin run of 10 m cell resolution MIN DEM (Fig. 10) and 5 m cell resolution TIN DEM by the used Sink Depth equals Sink Depth (max) values (Fig. 9). This could be taken as an indication of correct basin delineation results revealed by the natural mire surface pattern. From all Basin runs it could be concluded that, in spite of considerable variation between the basin delineation, and counts from different modelling runs, there were some general, prevailing delineation, which could be easily aggregated from Basin runs with higher basin count results. Therefore, the selection of one or other Basin run result, as the final decision option, depends on the set up of the scale of the study area extension, e.g. the whole mire landscape or specific ecotopes or microtopes.
In summary, we can draw the following conclusions: be used for preliminary mire surface study before ordering the LiDAR scanning for a non-studied mire landscape. 7. Scale-based accuracy and resolution of DEMs are a crucial issue in relation to simulation and visualization of the morphological and hydrological properties of the mire surface. 8. Since the detailed man-made topological and ecohydrological measurements of different mire complexes are usually labour-consuming, expensive and limited by the man-capability and skills, the results of remote sensing products must be used in mire landscape research. 9. The laser altimetry can be considered as the most accurate technology for the generation of detailed 3D DEMs for relatively flat but patterned mire surfaces.