EFFECT OF DIGITAL ELEVATION MODEL MESH SIZE ON GEOMORPHIC INDICES : A CASE STUDY OF THE IVAÍ RIVER WATERSHED – STATE OF PARANÁ , BRAZIL

Geomorphometry is the science of quantitative description of land surface morphology by the mean of geomorphic indices extracted from Digital Elevation Models (DEMs). The analysis of these indices is the first and most common procedure performed in several geoscience-related subjects. This study aims to assess the impact of mesh size degradation on different local and regional geomorphic indices extracted for GDEM and TOPODATA DEMs. Thus, these DEMs, having a mesh size of 30 m, were subsampled to 60, 120 and 240 m and then geomorphic indices were calculated using the full resolution DEM and the subsampled ones. Depending on their behavior, these indices are then classified into stable and unstable. The results show that the most affected indices are slope and hydrographic indices such as Strahler order, stream sinuosity and fractal dimension and watershed perimeter, whereas elevation remains stable. It also shows that the effect depends on the presence of the canopy and geological structures in the studied area.


Introduction
The analysis of the morphologic aspects of an area is one of the first and most common procedures carried out in hydrological, geomorphological, environmental and territorial planning studies.This analysis can be done through the calculation of several geomorphic indices that can be then used to identify morphologically homogeneous areas.These indices were used in many studies as an auxiliary tool such as in Sodré et al. (2007), Cherem et al. (2008), Bispo et al. (2009), Thommeret (2012) and Valeriano and Rossetti (2017).
The dissemination of global-coverage, free-access Digital Elevation Models (DEMs) supported several research areas related to geomorphological mapping, specifically the quantitative morphometric description and characterization, which has benefited greatly from the high resolution of these DEMs.The use of these DEMs for the derivation of morphologic descriptors have been used in different geomorphological analysis studies to infer interpretations and numeric analyses (Valeriano andRossetti 2011, 2012).
According to Pike, Evans and Hengl (2008), despite the fact that both classic morphometrybefore DEMand modern morphometry focus on spatial analysis, indices extraction and terrain surface analysis, a clear difference exists between both areas.Indeed, while classic morphometry especially addresses landform analysis, based on hypsometry, slope, frequency of elevation, elevation range, drainage density and topology, and terrain classification, modern morphometry focuses on improving elevation data processing, topography description and visualization, as well as numerical analysis based on digital data expressed by DEMs.
The quality of the geomorphic indices depends on the production technique and its parameters (El Hage et al. 2012), on the terrain morphology and its landcover (El Hage et al. 2017) as well as on the mesh size (Thompson, Bell and Butler 2001, Tang 2003, Kienzle 2004, Erskine et al. 2007, Vaze and Teng 2007, Wu, Li and Huang 2007, El Hage et al. 2010, Vaze, Teng and Spencer 2010, Polidori et al. 2012).The mesh size constitutes one of the main DEM characteristics and it is part of the triple scale concept manifested by dimension, extension and spacing, directly controlling the size of objects that could be extracted from DEMs (Blöschl and Sivapalan 1995).Given the importance of the scale on objects morphology representation, it is essential to study the impact of the mesh size on the information extracted for DEMs (Pike 2000, Pike, Evans andHengl 2008).The purpose of this study is to assess the sensitivity of the geomorphic indices extracted from DEMs to the change of the mesh size.Two DEMs are used, GDEM and TOPODATA, which are produced using two different techniques.First, the study area is presented.Then, the adopted datasets and methodologies are explained.Finally, the results are exposed and discussed.

Study area
The Ivaí River watershed (Figure 1 Due to its heterogeneous geologic, topographic and hydrologic characteristics, Destefani (2005) divided the Ivaí River into three sectors: upper, middle and lower (Figure 1).This subdivision was adopted for this study, and it was extended to the corresponding catchment areas.Thus, 58 subwatersheds, distributed along the three sectors and which have drainage areas greater than 100 km², were assessed.The upper sector is generally characterized by its convex and concave slopes, flat and encased peaks, narrow and deep valleys, and a higher degree of terrain dissection.The middle sector relies on predominant convex slopes, elongated peaks, encased valleys and a lower degree of dissection compared to the upper sector.The lower sector presents convex slopes, flat peaks, open valleys and low-degree dissection, and it is where the Ivaí river floodplain is located (Maack 2002, Santos et al. 2006).
The studied watershed presents a dendritic drainage pattern.In its upper and middle sectors, however, this pattern is conditioned by geologic structures, forming rectangular channels parallel to those structures (Figure 8) (Maack 2002).The Ivaí River watershed streams, as a whole, have a winding pattern with irregular contours and sudden breaks.According to Destefani (2005), the stream sinuosity shows both higher or lower amplitude and frequency, depending on the local geology and geomorphology.

Datasets
The Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) is one of five remote sensors onboard Terra Satellite, launched by NASA in 1999.The sensor carries three telescopic subsystems that enable the acquisition of images with different spatial resolutions depending on the wavelength range.For instance, the visible and near infrared (VNIR) subsystem has a 15 m pixel size, the short wave infrared (SWIR) has a 30 m pixel whereas the thermal infrared (TIR) has a 90 m pixel (Abrams 2000, Abrams andHook 2001).Each ASTER scene covers an area of 60 x 60 km (Abrams 2000).
The Global Digital Elevation Model (GDEM) is derived from the VNIR subsystem using nadir and backward images.It has a 30 m mesh size and it was created by processing 1.3 million ASTER scenes, covering the areas situated between 83°N and 83°S latitudes (Abrams andHook 2001, ERSDAC 2005).These scenes are used to generate relative DEMs (without ground control points) having a vertical accuracy of 10 m, and absolute DEMs (with ground control points) having a vertical accuracy of 7 m (Abrams and Hook 2001).
TOPODATA is a DEM obtained by oversampling the 3 arc-second Shuttle Radar Topography Mission (SRTM) DEM over Brazil.The SRTM mission, which took place in the year 2000, aimed to collect altimetric data for 80% of the Earth's land surface using SAR interferometry operating at X (λ = 3,1 cm) and C (λ = 5,6 cm) bands.This mission resulted in a global DEM with a mesh size of 1 arc-second (~30 m) for the United States and 3 arc-seconds (~90 m) for the other areas, having a vertical accuracy of 16 m and a horizontal accuracy of 20 m with 90% confidence (Rabus et al. 2003, Ludwig and Schneider 2006, Berry, Garlick and Smith 2007).TOPODATA is generated by oversampling SRTM-DEM from 90 m to 30 m using an algorithm based on kriging, a geostatistical interpolation method.The extraction methodology is further explained in Valeriano and Albuquerque (2010) and Valeriano and Rossetti (2012).
In addition to these DEMs, 1:50,000 and 1:100,000 topographic maps are used.These maps are generated basing of photogrammetric surveys conducted by IBGE (Brazilian Institute of Geography and Statistics) and DSG (the Ministry of the Army -Directorate of Geographic Services) between 1960 and 1970.The maps were then provided free of charge by the ITCG (Institute of Land, Cartography and Geosciences) for the Brazilian southern state of Paraná.

Methodology
To assess the impact of the mesh size on the geomorphic indices, the DEMs should first be subsampled.Thus, GDEM and TOPODATA are subsampled by keeping one point every two points.Therefore, the first subsampled DEM will have a mesh size of 60 m.Then, the subsampled DEM in his turn is resubsampled by the same method (one point every two), so that its new mesh size is 120 m and 240 m for the last one.This procedure is conducted with the resampling tool of ArcGIS 10.1, and the nearest neighbor interpolation is used, which assign the elevation of the closest pixel to the output pixel.Then, hydrographic and geomorphic indices are extracted from the full resolution DEMs and the subsampled ones and from the topographic maps in order to compare them.
The drainage network extraction and watershed delineation are performed using ArcGIS 10.1, in which the flow direction is determined using D8 algorithm and the extraction chain involves the following steps: 1. First, the DEM is filled where sinks and peaks that block the flow are removed.Such procedure avoids possible mistakes in the generated drainage network.
2. Then, the flow direction and flow accumulation are extracted.
3. Finally, the stream order is calculated and the sub-watersheds are delineated.
Initially, the geomorphic indices were classified into local and regional.Local indices characterize the value at a specific point, such as elevation and slope, whereas regional indices are descriptors of an area of the watershed, such as the form factor (Table 1).Form factor or conformation index regional Ratio of the basin area to the square of the basin length The calculations have been carried out by ArcGIS and SSM (Self-Similarity Map) software (developed by INRA (Institut National de la Recherche Agronomique)).The latter is dedicated to quantify the geometry and topology of stream networks and, therefore, it enables the quantitative analysis of any channel network (Gaucherel, Salomon and Labonne 2011).

Results
Elevation histograms (Figure 2) of the studied watershed show a higher number of pixels for the same elevation value for GDEM compared with TOPODATA.These small differences are related to the production techniques and methods of these DEMs (Polidori, El Hage and Valeriano 2014).This figure also shows that the increase of the mesh size does not affect the elevation histogram for the both DEMs.The slope histograms (Figure 3), show that an increase in the mesh size induces a disappearance of steep slopes and low slope angles become more abundant, which is an expected result since the terrain surface becomes smoother.This can have an important effect in studies that search to map unstable areas.Indeed, the extent of unstable regions that are characterized by the presence of steep slopes, tend to be underestimated with the increase of the mesh size as these slopes disappear (El Hage et al. 2010, Polidori el al. 2012).Moreover, the impact of subsampling is more pronounced on GDEM than on TOPODATA.This could be due the strong effect induced by the presence of canopy in GDEM that increases the number of steep slopes in the DEM, specifically with small mesh sizes (Polidori and Simonetto 2014).This effect diminishes with the increase of the subsampling ratio, and with a mesh size of 240 m, the histograms of the two DEMs become very similar.To assess the impact of mesh size on the hydrography, stream networks were extracted from both DEMs.Some studies evaluated the drainage network extracted from GDEM and they concluded some imperfections in the stream morphometry (Hosseinzadeh 2011) and a good positional accuracy of the streams (Polidori et al. 2014).Figure 4 shows that there is a loss in the topologic and morphologic details of streams obtained from GDEM and TOPODATA compared with the drainage extracted from the topographic map at a scale of 1: 50,000.The difference is more significant in the upper and lower sectors for both GDEM and TOPODATA.
Likewise, the subsampling has also an effect on the drainage network for both DEMs, and this is visually clear as shown in Figure 5.This impact is manifested by a reduced number of streams having small orders and a decrease in the morphological details with slight differences between GDEM and TOPODATA.Figure 6 shows the relation between Strahler order (Strahler 1957) and the number of streams in the lower sector for the different mesh sizes.The streams extracted from full resolution GDEM are more numerous the those extracted from full resolution TOPODATA for all orders, which also could be an effect induced by the canopy.However, as the mesh size increase, the number of streams decreases for all the orders and for both DEMs, and the difference between the two DEMs becomes smaller.In general, the sinuosity of the streams in Ivaí River watershed is moderate, and their extent and frequency depend on local geology and geomorphology.Figure 7 represents the mean value of sinuosity and fractal dimension for the different Strahler orders and different mesh sizes.This figure shows a big difference between the different sectors with respect to the mesh size.Basing on Figure 8, a relation can be inferred between this difference and the presence of a complex system of fractures in the basin.Indeed, the upper sectors is characterized by a very dense system of structural lineaments with three main directions: NE, NW and EW (the lineaments are extracted using a method based on satellite imagery and geophysical measurements), the middle sector has a lower density and lower sector has the least density of lineaments that do not have a main direction (Wildner et al. 2006).
The sinuosity and fractal dimension of the first two Strahler orders increases with the increases of mesh size and decrease for the higher orders in the upper and middle sectors.However, in the lower sector, the increase of mesh size does not significantly affect these indices (Figure 7).The low value of sinuosity and fractal dimension for the first two orders in the case of full resolution DEMs is due to the presence of an intensive network of faults in the upstream of the upper and middle sectors (Figure 8).Indeed, these faults are linear objects and their sinuosity and fractal dimension values tend to be 1.An increase in the mesh size causes deformation of these linear structures and then an increase in their sinuosity and fractal dimension.The confluence angles in the upper and middle sectors are conditioned by local geology as shown in the Figure 9, where it is possible to visualize the influence of lineaments in the streams.Figure 10 represents the confluence angles generated with different mesh sizes.In general, the mean value of confluence angles decreases when the mesh size increases.This is partly due to the intensive presence of linear structures that with NE-NW direction that create confluence with near right angles (this is why the mean value is high).The subsampling produces a smoothing in the relief and the angles between these structures tend to decrease as a result.The drainage area shows slight difference between full resolution and subsampled DEMs and also when compared to topographic maps (Figure 11 -A).A small difference is shown in the lower sector, where the terrain is very rugged and the subsampling can cause the disappearance of some terrain details and consequently the sub-basin boundaries differs (Figure 12).
By comparing the values of the perimeters extracted from topographic maps and DEMs for all the watershed sectors studied (Figure 13), an inverse relation that exist between the mesh size and the perimeter.Indeed, the more the mesh size, the less the perimeter and this is due to a decrease in the surface roughness, where noise and small landforms disappear.It is important to highlight that the perimeter values for the 30 m mesh size in both DEMs are slightly higher than those of the topographic map.This may be due to the presence of canopy that has a clearer effect on 30 m DEM than on the subsampled ones (Polidori and Simonetto 2014), specifically in the case of GDEM.Moreover, the perimeters calculated by GDEM are always higher than those calculated by TOPODATA except the last mesh size, which also implies the effect of the canopy that tends to disappear with a mesh size higher than 120 m.The form factor (Figure 11 -B) undergoes a small change in value with the increase of the mesh size.For example, the drainage area in the upper sector is almost stable for all DEMs, but the form factor presents a difference in value, which comes from a difference in the perimeter of the subwatershed (Figure 13).

Discussion and conclusion
This study aimed to assess the impact of subsampling on different geomorphic indices extracted from DEMs.Local and regional indices were evaluated and the results show that the mesh size produces a significant effect on some indices like slope and watershed perimeter and low effect on other indices like drainage area and elevation.Table 2 represents the obtained stable and unstable indices.
The effect of mesh size depends on presence of canopy in the studied area.Indeed, the canopy modify the morphology on the DEM and the geomorphic indices that arise from are affected specifically with a small mesh.According to Polidori, El Hage and Valeriano (2014) and that is why the effect of the canopy is low in TOPODATA.The increase of the mesh size decreases the effect of the canopy in GDEM and with a mesh size higher than 120 m, the difference between GDEM and TOPODATA becomes very small.The impact also depends on the presence of structures in the studied area.Indeed, the lineaments are better represented in full resolution DEMs, where a subsampling tends to modify these linear features, which affects the geomorphic indices, specifically the morphometry of the stream network.Strahler order (regional) Form Factor (regional) Sinuosity (regional) Fractal dimension (regional) Perimeter (regional) The drainage network extracted from the 30 m GDEM and TOPODATA showed a difference in the stream network compared with that extracted from the topographic map.In the floodplain areas, the movement of water is hardly detected using a DEM as it is controlled by very small topographic features, specifically if the mesh size is not fine enough to represent these small features such as narrow channels connecting main rivers and floodplains (Yamazaki et al 2012, Saksena andMerwade 2015).Moreover, the presence of noise and canopy in the DEM highly affect the drainage network in these regions, specifically with small mesh sizes (Lane, James andCrowell 2000, Polidori andSimonetto 2014).
According to Fernández (2011), the studied DEMs have limitations in detailed drainage network extraction (scale > 1:50,000).With TOPODATA, though, it is possible to generate coherent drainage network compatible with a scale of 1:100,000; although it still requires editing work, especially in flat areas as floodplains because of a limitation in the drainage extraction algorithms.Nevertheless, in low slope angles, GDEM and TOPODATA are compatible with scale of 1:50,000, whereas in steeper slopes, like the upper and middle sectors, the compatible scale is 1:100,000 or smaller (Moura, Bias and 2014).
Flow direction extraction algorithm (D8) has some limitations and may, therefore, cause some errors.Algorithms such as D8 (unidirectional) gives better results on converging topographies.In diverging topographies, however, it generates errors since it restricts the direction of the flow towards a single pixel among its eight neighboring pixels (Fernández 2011, Fernández et al. 2012).
There is also an influence caused by other factors such as DEM elevation error that propagates in the geomorphic indices (Tang et al. 2001, El Hage et al. 2017).
It is important to emphasize that the watershed in question is strongly controlled by a series of structural lineaments, which is prominently reflected on the hydrographic network and, consequently, on the local and regional indices, especially in the upper and middle sectors (Santos 2015).Another watershed located in a smoother terrain and missing intense structural lineaments could provide different results.Further research is needed to profoundly understand the impact of mesh size on the geomorphic indices over these structures.
In this study, the influence of scale change on geomorphic indices has been assessed.However, to have a better idea on how the mesh size change affect the geomorphology of the DEM, the quality of the aforementioned indices should be assessed using ground control data, and this will be investigated in future research.It is also important to point out that the results obtained depend on the extraction algorithms adopted in the study, where other algorithms would give slight differences.
), located in the State of Paraná (South Brazil), is the second largest watershed in the state, consisting of an area of about 36,533 km 2 , which corresponds to 18% of Paraná territory.The Ivaí River is an important left bank tributary of the Paraná River formed by the confluence of Patos and São João rivers (located on Serra da Boa Esperança, municipality of Ivaí), at an elevation of about 480 m, and it flows at a distance of about 800 km towards its estuary in the Paraná River (Pontal do Tigre, municipality of Icaraíma), at an elevation of about 240 m (Destefani 2005, Meurer, Bravard, and Stevaux 2011).

Figure 1 :
Figure 1: Location of the Ivaí River watershed, its topography and the three sectors.

Figure 2 :
Figure 2: Elevation histograms of the full resolution GDEM (A) and TOPODATA (B) DEMs and the subsampled ones.

Figure 3 :
Figure 3: Slope histograms of the full resolution GDEM (A) and TOPODATA (B) DEMs and the subsampled ones.

Figure 4 :
Figure 4: Drainage network extracted from topographic maps and from GDEM and TOPODATA in (A) the upper sector; (B) the middle sector (C) and the lower sector.

Figure 5 :
Figure 5: Drainage networks extracted from GDEM (A) and TOPODATA (B) in the lower sector.

Figure 6 :
Figure 6: Number of streams for each Strahler order in the lower sector.

Figure 7 :
Figure 7: Mean sinuosity and mean fractal dimension for all Strahler orders in the upper sector (A) the middle sector (B) and the lower sector (C).

Figure 9 :
Figure 9: River confluences in Ivaí River watershed in the upper sector (A) the middle sector (B) the lower sector (C) and the floodplain area (D).

Figure 10 :
Figure 10: Mean values of confluence angles extracted from GDEM and TOPODATA.

Figure 11 :
Figure 11: Drainage area (A) and form factor (B) calculated for topographic maps, GDEM and TOPODATA DEMs for the upper, middle and lower sector sub-watersheds.

Figure 12 :
Figure 12: Boundaries of a sub-basin extracted from GDEM (A) and TOPODATA (B) in the lower sector.

Figure 13 :
Figure 13: Perimeter of the sub-basins extracted from the topographic map, GDEM and TOPODATA in (A) upper sector; (B) middle sector (C) lower sector of the watershed.

Table 1 :
Local and regional geomorphic indices used in the study.
, TOPODATA carries the advantage of extending SRTM geomorphologic properties to a rather local scale, since the oversampling method of SRTM from 90 m to 30 m is based on geomorphologic assumptions Effect of ...

Table 2 :
Stable and unstable geomorphic indices.