Influence of cell size on volume calculation using digital terrain models: A case of coastal dune fields

In this work, we analyze how variation in cell size influences the volume calculated from Digital Terrain Models (DTMs) derived from a LiDAR (Light Detection and Ranging) survey in two coastal Late Holocene dune fields in southern Brazil. Cell size varied from 1 to 100 m. The RMSE (Root Mean Square Error) of the resampled DTMs from the original LiDAR (with 0.5 m resolution) increases linearly with cell size, while the R-squared decreases following a secondorder trend. The volume does not show simple linear or exponential behavior, but fluctuates with positive and negative deviations from the original DTM. This fluctuation can be explained by a random factor in the position of the cell with regard to landforms and a relationship between cell and landform size, wherein a small change in cell size can lead to an underor overestimation of volume. ASTER GDEM (Global Digital Elevation Model) and X-SAR SRTM (Shuttle Radar Topography Mission) 1 arcsec Digital Elevation Models (DEMs) were not considered viable volume sources due to large deviations from the reference data, either as a consequence of noise in the SRTM X-SAR data or lack of bias elevation correction to a common reference base in the GDEM processing chain. Volumes from a 3-arcsec SIR-C SRTM deviated around ± 5% from the reference data and are considered suitable input for numerical simulations of Quaternary dune field evolution models because these values should be within the expected range of sediment volume changes over hundreds to millions of years.


Introduction
Volume calculation from Digital Elevation Models (DEMs) or Digital Terrain Models (DTMs) 1 is a very important tool in the study of Quaternary landforms. Aeolian dune fields stand out as widespread landforms occurring in wet to dry depositional settings (Fryberger and Dean, 1979;Short, 1988;Wang et al., 2002;Livingstone et al., 2007;Martinho et al., 2010) on Earth as well as on other planetary bodies such as Mars, Venus and Saturns moon Titan (Hayward et al., 2007;Radebaugh et al., 2008;Bourke et al., 2010;Fenton and Hayward, 2010). The sand supply for dune fields formation and growth is sensitive to shifts in wind patterns and pluviosity. Thus, changes in dune field volume and morphology can be related to climate change (Clemmensen et al., 2007;Sawakuchi et al., 2008;Tsoar et al., 2009;Singhvi et al., 2010;Levin, 2011). Despite the surface morphology of a dune field provide important information about wind patterns and sand availability, it records only the final phase of the evolution of a dune field. In a source-to-sink system approach, dune fields are reservoirs of sand supplied by other depositional systems such as a wave-dominated coastal system or a fluvial system. The volume of sand stored in a dune field depends on the sand flux through time, which is sensitive to climate change. Thus, the calculation of dune field volume combined with dating methods to determine ages of aeolian sand deposition, such as that based on cosmogenic nuclides (Vermeesch et al., 2010) or luminescence (Wintle, 1993;Lancaster, 2008), permit to obtain a more complete record of the climate variables affecting the evolution of dune fields and their coupled depositional systems. Numerical simulations play an important role in understanding the response of dune fields to climate change. In this context, numerical simulations of dune field evolution can be greatly improved if the volume of sediment within a certain area is known within an acceptable deviation from its true value.
In this work, we analyze how the spatial resolution of DTMs derived from a high-resolution LiDAR (Light Detection and Ranging) dataset influences the output of volume calculations. Given the high cost of LiDAR surveys and the large areas occupied by dune fields, it is also interesting to evaluate whether coarser-resolution DEMs, such as SRTM (Shuttle Radar Topography Mission - Farr et al., 2007) or ASTER-GDEM (ASTER Global Digital Elevation Model - Abrams et al., 2010), can be used to estimate the volume of dune fields within an acceptable deviation range. The calculation of sediment volumes using SRTM DEMs would be very useful for studying Quaternary aeolian sedimentation in Brazil, where more than two hundred dune fields are recognized (Sawakuchi, 2006;Giannini et al., 2007), as well as on a global scale.

Study Area
The study area is located in Santa Catarina State, southern Brazil ( Fig. 1). This area comprises barrier-lagoon depositional systems with associated dune fields (Giannini et al., 2007). These systems evolved since the Middle Holocene during a phase of low-rate decreasing relative sea level after a highstand of around 2.5 m above the present sea level (Angulo et al., 2006). The Garopaba (northern sub-area) and the Itapirubá (southern sub-area) dune fields comprise unvegetated and vegetated aeolian dunes. The unvegetated dunes are represented by mostly barchanoid chains, while the vegetated ones include parabolic dunes, blowouts and foredunes (Martinho et al., 2006;Giannini et al., 2007;Hesp et al., 2007Hesp et al., , 2009). These dune fields are a result of wind strength intensification and sand supply increase in southern Brazilian coast during the Late Holocene. The active portion of both dune fields, which are formed by undisturbed barchanoid chains, present ideal conditions for studying the relationship between cell size and volume estimation on DEMs.

LiDAR processing
LiDAR data were collected on October 2010 by Geoid Laser Mapping using an Optech ALTM 3100 sensor with a saw-tooth scanning pattern, a density of about one point per 0.5 m, measured from an altitude of 1200 m (4000 ft). Raw LiDAR data (with up to four laser pulses) were processed by Geoid using proprietary algorithms, and bare ground data (i.e., filtered from vegetation and man-made structures) were delivered with a vertical accuracy of 0.15 m (1σ) and a horizontal accuracy of 0.5 m (1σ).
Ground point data were converted to ASCII xyz points with ASPRS libLAS 1.6 (Butler et al., 2011) and imported into GRASS-GIS (Neteler and Mitasova, 2008;GRASS Development Team, 2009) as vector points. A DTM representing the ground surface with 0.5 m spatial resolution was created by interpolation of the vector points with bicubic splines .

SRTM and ASTER GDEM data
To evaluate whether freely available coarser-resolution DEMs could be used to estimate the volume of dune fields within an acceptable deviation range, 3-arcsec SRTM SIR-C (Farr et al., 2007), 1-arcsec SRTM X-SAR (Eineder et al., 2001) and 1-arcsec ASTER GDEM (Earth Remote Sensing Data Analysis Center, 2009; Abrams et al., 2010) data were also imported into GRASS. Volumes were calculated with a pixel size of 90 m for 3-arcsec data and 30 m for 1-arcsec data.

DTM resampling and volume calculation
To calculate the sand volume, we used the GRASS-GIS r.volume module (Hinthorne, 1988). This module calculates volume by summing cell values within a given area and then multiplying by the area occupied by those cells. An elevation of 0 m (zero) was a reference base level. Although a 3 horizontal surface does not represent the true boundary of the dune field bottom surface, which can be irregular due to underlying paleotopography (Giannini et al., 2007), this should not interfere with the goals of this paper because we are interested in a comparison between volume values at different spatial resolutions. The bottom surface of aeolian dunes can be detected with geophysical methods, such as Ground Penetrating Radar (GPR). The Garopaba and Itapirubá dune fields migrate over mostly low relief (< 2 2.5 m) beach ridges. The similarity between the sediments of beach ridges and aeolian dunes makes it difficult to distinguish between the two with GPR, rendering this approach impractical. One point to consider is a DEM that has areas with positive and negative values. A negative elevation value results in a volume with negative value (sign), which will be subtracted from the positive volume when the whole DEM is processed. In this work, we are interested in sand volume deposited over the coastal area above mean sea level; therefore, we discarded any negative elevations before calculating volumes.
The 'true' volume was taken as the volume of the DTM with a 0.5 m resolution. To evaluate the effect of spatial resolution over volume values, the cell size varied from 1 to 100 m, and at each step, the original DTM was resampled and the volume was calculated. Resampling was achieved by calculating the mean elevation value for a given resolution. This resampling simulates data acquired by sensors with different spatial resolutions; therefore, it is appropriate that the elevation value is the average of the actual elevations within the cell (Grohmann et al., 2010).
For each resolution step, using all cell values of the resampled DTM, we calculated the volume difference (deviation), the Root Mean Square Error (RMSE, absolute error) and the Pearson correlation coefficient (R-squared, goodness of fit) between the resampled DTM and the LiDAR DTM , where and where z is the elevation for each cell and n is the total number of cells. Additionally, topographic profiles were derived from all of the models to visualize how the variation in spatial resolution would affect the representation of landforms. Fig. 2 shows the calculated volume of the Garopaba and Itapirubá dune fields for each resolution step and the volume difference between the models as a percentage of error from the LiDAR DTM.

LiDAR DTM
The maximum deviation in volume ranges from ±5% for the Garopaba field and ±4% for the Itapirubá field, with a significant high frequency alternation between positive and negative values for cell sizes above 20 m. This result is somewhat counterintuitive, as one could expect the difference in volume to increase with cell sizes. In fact, the RSME and R-squared between the resampled models and the LiDAR DTM present such behavior with an approximate linear increase. In Fig.  3, the left plots show a linear trend of increasing RSME values as the cell size increases, and the right plots show a decrease in R-squared values that can be fit to a second-order polynomial trend.
The linear relationship between RSME values and cell size shown in Fig. 3 should be seen as an approximation for cell sizes less than 10 m, as the data points show non-linear behavior in this area of the plot, particularly for cell sizes under 5 m. This is an important point for future studies, as LiDAR technology keeps evolving toward finer resolutions. 4  From these profiles, one can see that the final volume can be over-or underestimated, depending on the size of the cell, the size of the landforms and the position of the cell with reference to the landforms. Considering the resampled DTMs, the profile for 40 m cells shows a fair agreement with the reference surface and its calculated volume (10.67 × 10 6 m 3 ) deviates 0.60% from the true LiDAR DTM (10.60 × 10 6 m 3 ). A cell size of 80 m gives an underestimation of the volume (10.42×10 6 m 3 , -1.70 % from the reference), while a resolution of 87 m is sufficient to overestimate the volume (10.97 × 10 6 m 3 , +3.40 % from the reference). 5

SRTM and ASTER GDEM
Because the study area comprises active dune fields and there is a 10 year difference between SRTM and LiDAR data, the migration rate of these landforms must be considered before comparing calculated volumes. The SRTM Mission was flown on February 2000 (Farr et al., 2007), while ASTER-GDEM was created from optical imagery collected between 2000 and 2007 (Abrams et al., 2010) and the LiDAR data were surveyed in October 2010.
A migration rate of approximately 3 m/year toward the S-SW is reported by Giannini et al. (2005) for an area immediately north of Garopaba. The time span between SRTM and LiDAR datasets would account for a 30 m shift in the dune's position, that is, one pixel in SRTM X-SAR and ASTER GDEM and 1/3 pixel in SRTM SIR-C. We do not consider such a difference as an important factor in our analysis. Table 1 shows the volumes calculated from the LiDAR DTM, SRTM SIR-C, SRTM X-SAR and ASTER GDEM, as well as the differences between these products and the reference LiDAR data.
The deviation of SRTM SIR-C data (90 m) from the LiDAR DTM was of -5.06% for the Garopaba dune field and of +5.46% for the Itapirubá dune field. SRTM X-SAR (30 m) presented a higher variability, deviating +9.79% at Garopaba and +4.08% at Itapirubá, while ASTER GDEM (30 m) had the highest deviation range, from +36.59% at Garopaba and only +0.32% at Itapirubá.
Topographic profiles for the Garopaba (Fig. 5) and Itapirubá (Fig. 6) dune fields provide a visual comparison between the original LiDAR DTM, a DTM resampled to 90 m and SRTM SIR-C (Figs. 5A,6A), as well as between SRTM X-SAR and ASTER GDEM (Figs. 5B,6B). 6  Fig. 5A, in the Southwest area of the Garopaba dune field (left side of plot), SRTM SIR-C overestimates the landforms, while in the Northeast (right side), where the altitude is very low, it underestimates them, likely because of interference with the radar signal due to the presence of water. In Fig. 5B, the SRTM X-SAR shows a strong noise component, with altitude values ranging from +60 m to -22 m and a poor representation of the topography. The ASTER GDEM data proved to be consistently above the reference topography. This could be due to the lack of bias elevation correction to a common reference base in the GDEM processing chain (Reuter et al., 2009). Profiles for the Itapirubá dune field show that SRTM SIR-C (Fig. 6A) presents a general representation of the topography, with a tendency of overestimation of elevation, more clear in the central area. In this case, the ASTER GDEM (Fig. 6B) provides a good picture of the topography, but the SRTM X-SAR still has a strong noise component.

Discussions and Conclusions
In this work, the analysis of volume variation with cell size showed a fluctuation of calculated volume values from the reference LiDAR DTM, instead of a simple increasing or decreasing behavior. This can be explained by a random factor between the size of the cell, the size of the landforms and the position of the cell with reference to the landforms, where a small change in cell size can lead to an under-or overestimation of the volume.
An optimal cell size of 20 m for the Garopaba dune field could be interpreted from the plots in Fig.2 because this value represents a break point, after which deviation from the reference DTM becomes larger and starts to fluctuate. From the plots of Fig.3, however, one could interpret an optimal size of 40 m, after which both RMSE and R-squared for the LiDAR DTM drop more rapidly. The volumes for the Itapirubá dune field show steady behavior up to 25-30 m (Fig.2B), but the RMSE and R-squared for the LiDAR DTM start dropping only after 50 m.
The size of the landforms must also be considered. In the Garopaba dune field, despite an altitude of approximately 40 m at its Southwest area (left side of plots in Fig. 5), which could be related to the paleotopography, the height of dunes remains under 10 m. In the Itapirubá dune field, dunes exceed 20 m in height (Fig. 6).
Although the volume of dune fields is an important variable in numerical simulations to study dune field dynamics, it is hard to determine in a fast and inexpensive way considering the existence of hundreds of dune fields along the Brazilian coast. The high cost of LiDAR surveys is still a barrier to this type of data, especially if one intends to study large areas such as desert dune fields. 8 Numerical models are necessary to study the evolution of aeolian dune fields whose development occurs over a timescale of hundreds to millions of years (Livingstone et al., 2010;Vermeesch et al., 2010) because this timescale is beyond observational data. Numerical simulations address variations of the dune field properties (volume, area, etc.) over time spans from days to hundreds or million of years. In such scenarios, a fluctuation in sand volume of 5.5% would be greater than that expected for short periods (days to years), but not greater than that expected for longer periods of time (hundreds to millions of years) (Sawakuchi, 2006).
With a deviation of volume from the LiDAR DTM of around +/-5.5% for both dune fields, SRTM SIR-C (90 m resolution) is a viable source for volume values of dune fields. We also expect this deviation to decrease as the size of dunes and dune field volume increases, as suggested by the comparison between the Garopaba and Itapirubá dune fields.
The ASTER GDEM still generates some controversy in regard to its usefulness. Although it has been found useful in the study of large aeolian landforms (Bubenzer and Bolten, 2008;Hayakawa et al., 2008;Hugenholtz and Barchyn, 2010), Reuter et al. (2009) note that the GDEM product should not be used in terrain analysis at a 30 m resolution, and Miliaresis and Paraschou (2011) stress the importance of the number of scenes (stack number) used to derive the DEMs, and how vertical accuracy could be elevation-and slope-dependent. In our case, ASTER GDEM data for the Garopaba dune field proved to be consistently above the reference topography, resulting in a deviation of more than 35% from the reference data, which could be related to unresolved issues in the GDEM processing chain. We also note that the Garopaba and Itapirubá dune fields are relatively small landforms if compared to huge desert fields, such as the Namib dune field where the ASTER data proved to be reliable estimates of sand volume (Bullard et al., 2011).
The SRTM X-SAR data presented a strong noise component and very poor representation of the landforms in the Garopaba dune field. According to Eineder et al. (2001), noise in the 30 m scale for the X-SAR SRTM is due to thermal effects and should be in the ± 4 m range, but in Fig  5, the X-SAR elevation ranges from +60 m to -22 m. We believe these issues could be related to radar shadowing in the steep Southwest-facing slopes of the dunes (Strozzi et al., 2003;Ludwig and Schneider, 2006).
The strong difference in the behavior of the SRTM X-SAR between the dune fields, with large deviations from the LiDAR DTM at Garopaba and only small differences at Itapirubá, could be attributed to the geometry of the dune fields, where the Garopaba landforms favor radar shadowing effects, while the errors of the ASTER GDEM show the need for attention when working with such data. A simple comparison to 3-arcsec SRTM data using topographic profiles is recommended prior to analysis.
The high periodicity and geometry of these aeolian landforms may have influenced the results we obtained. Other more random landscape configurations may present a different behavior, and we expect to test such a hypothesis in future studies.