Hydrol. Earth Syst. Sci., 14, 393–405, 2010 www.hydrol-earth-syst-sci.net/14/393/2010/ © Author(s) 2010. This work is distributed under the Creative Commons Attribution 3.0 License. Hydrology and Earth System Sciences Reading the bed morphology of a mountain stream: a geomorphometric study on high-resolution topographic data S. Trevisani, M. Cavalli, and L. Marchi CNR IRPI, Corso Stati Uniti 4, 35127 Padova, Italy Received: 11 November 2009 – Published in Hydrol. Earth Syst. Sci. Discuss.: 27 November 2009 Revised: 4 February 2010 – Accepted: 22 February 2010 – Published: 26 February 2010 Abstract. High-resolution topographic data expand the potential of quantitative analysis of the earth surface, improving the interpretation of geomorphic processes. In particular, the morphologies of the channel beds of mountain streams, which are characterised by strong spatial variability, can be analysed much more effectively with this type of data. In this study, we analysed the aerial LiDAR topographic data of a headwater stream, the Rio Cordon (watershed area: 5 km2 ), located in the Dolomites (north-eastern Italy). The morphology of the channel bed of Rio Cordon is characterised by alternating step pools, cascades, and rapids with steps. We analysed the streambed morphology by means of ad hoc developed morphometric indices, capable of highlighting morphological features at a high level of spatial resolution. To perform the analysis and the data interpolation, we carried out a channel-oriented coordinate transformation. In the new coordinate system, the calculation of morphometric indices in directions along and transverse to the flow direction is straightforward. Three geomorphometric indices were developed and applied as follows: a slope index computed on the whole width of the channel bed, directional variograms computed along the flow direction and perpendicular to it, and local anomalies, calculated as the difference between directional variograms at different spatial scales. Directional variograms in the flow direction and local anomalies have proven to be effective at recognising morphologic units, such as steps, pools and clusters of large boulders. At the spatial scale of channel reaches, these indices have demonstrated a satisfactory capability to outline patterns associated with boulder cascades and rapids with steps, whereas they did not clearly differentiate between morphologies with less marked morphological differences, such as step pools and cascades. Correspondence to: S. Trevisani (sebastiano.trevisani@irpi.cnr.it) 1 Introduction Channel bed morphology, resulting from the interactions between hydrology, topography and sediment supply, is sensitive to climate changes, land use changes, and to disturbances (natural or anthropogenic) in the upstream river network. Accordingly, the study of channel morphology provides information on the influencing processes mentioned above. Small high-gradient streams in mountainous areas have a complex and dynamic morphology, strongly influenced by the topographic ruggedness of the basin and the spatial variability of sediment supply. Geomorphic processes, such as debris flows and landslides are widespread and frequent in mountainous basins, and they often interact with the channels. These extrinsic episodic processes have a quite immediate effect on local bed morphology. The characterisation of topography and morphology of mountainous channels, the assessment of changes caused by large floods and the objective comparison between different channels require reliable data and procedures. However, quantitative studies of mountain streambed morphology are challenging both with regard to (1) the acquisition of extensive high-resolution topographic data and (2) difficulties in classifying complex stream bed morphology. 1. The topographic survey of mountainous streams requires collecting large amounts of data (Grant et al., 1990; Wohl et al., 1997; Chin, 1999; Lenzi, 2001; Milzow et al., 2006). Surveys carried out using traditional techniques, such as theodolite and stadia, may lead to high-quality representations of channel-bed topography, but are expensive and time-consuming. This can be overcome by means of Airborne LiDAR technology, which makes it possible to quickly collect large amounts of high-resolution topographic data, allowing careful surveys of large areas at reasonable costs and enabling the implementation of high-resolution DTMs. Published by Copernicus Publications on behalf of the European Geosciences Union. 394 S. Trevisani et al.: Reading the bed morphology of a mountain stream 2. The morphology of channel beds is characterised by a hierarchical structure with features that can be recognised at both the spatial scale of morphological units (e.g., steps and pools) and of stream reaches (e.g., step pool sequences, cascade, and plane bed). The classification of stream reaches into unequivocally defined classes is difficult when the channel bed morphology is complex, with short channel reaches of different morphologies and fuzzy transitions between them. A geomorphometric approach, i.e., the quantitative study of morphology (Pike, 2000; Hengl and Reuter, 2009), can provide valuable elements for interpreting the complex morphologies of mountainous streambeds and for handling the subjectivity intrinsic in the classification of bed morphologies by means of field observations. Airborne LiDAR data have been widely applied in many geomorphological fields, especially for studying landslides (Glenn et al., 2006; Ardizzone et al., 2007) and sediment dynamics (Scheidl et al., 2008; Brown et al., 2009). In studies on riverine environments, early applications of LiDAR data were mainly intended to provide improved DTMs for applying numerical flood models. Recent applications of aerial LiDAR data to the analysis of fluvial environments include the geomorphological mapping of gravel-bed rivers (Charlton et al., 2003), tidal channels (Mason et al., 2006) and river valley environments (Jones et al., 2007), the evaluation of river bank erosion (Thoma et al., 2005) and long-term changes caused by fluvial and debris-flow activity (Magirl et al., 2005), the assessment of the slopes in headwater channels (Vianello et al., 2009), the hillslope to valley transition (Tarolli and Dalla Fontana, 2009) and the delineation of the water surface (Höfle et al., 2009). In this study, morphometric indices derived from terrain analysis (Wilson and Gallant, 2000) and geostatistics (Goovaerts, 1997) are used to characterise channel bed morphologies through the analysis of a high-resolution DTM (cell size 0.5 m) developed from aerial LiDAR data. We read the streambed morphology using ad-hoc developed morphometric indices capable of highlighting morphological features at a high level of spatial resolution. Findings arising from the analysis of the morphometric indices are checked against field evidences to evaluate the correspondence with the morphological units classified by means of field observations. 2 Study area The study area is the main channel of the Rio Cordon, a small stream in the Dolomites (Eastern Italian Alps). The drainage basin of the Rio Cordon (Fig. 1) covers 5 km2 ; the elevation ranges between 1763 and 2748 m a.s.l., with an average value of 2200 m a.s.l.; the average slope is 27◦ . Slopes that locally exceed 45◦ are common, and subvertical cliffs are present in the upper part of the basin. The geological Hydrol. Earth Syst. Sci., 14, 393–405, 2010 Fig. 1. Geographical location and shaded relief map of the Rio Cordon. The rectangle outlines the studied portion of the channel. setting of the basin is rather complex: dolomite crops out in the upper part, whereas in the central and lower parts volcanic conglomerates, sandstones and calcareous-marly rocks crop out; moraines, scree deposits and landslide accumulations are also widespread. The basin presents a mean annual rainfall of approximately 1100 mm; precipitation occurs mainly as snowfall from October to April. Runoff is dominated by snowmelt in May and June; floods mostly occur in the summer and early autumn. Vegetation cover consists mainly of mountain grassland (60%) and widespread shrubs (15%), while forest stands composed by spruce and larch are found only in the lower part of the watershed and occupy 6% of the total area. Unvegetated areas (bare rock and scree) are common in the upper parts of the watershed. No relevant artificial structures, such as check dams or channel lining are present along the channels. In this study, we consider the lower part of the main channel reach, starting from a rocky gorge in the central part of the basin, and ending at the basin outlet (Fig. 1). The analysed channel reach has a total length of 1625 m and ranges in elevation from 2098 to 1818 m (ellipsoid height) (Fig. 2). The drainage area varies from 2.23 km2 at the upstream end of the studied channel reach to 5 km2 at its downstream end. Lenzi (2001) presented a classification of the channel bed morphology of the Rio Cordon and reported the geometric parameters of step-pool structures, measured using a total station. The classification of Lenzi (2001), also followed www.hydrol-earth-syst-sci.net/14/393/2010/ S. Trevisani et al.: Reading the bed morphology of a mountain stream 395 Fig. 2. Longitudinal profile and plan view of the studied channel reach, showing the classification of bed morphology. in Cavalli et al. (2008), identified three morphologies in the studied channel reach: a cascade reach at the outlet of the rocky gorge, two riffle-pool sequences, separated by a short “mixed” reach, and an alternation of step-pool and “mixed” reaches in the lower part of the stream. Although no relevant changes in the morphology have occurred in the Rio Cordon after the surveys reported by Lenzi (2001), the analysis of bed morphology was revised to achieve a more detailed classification of channel reaches. The revised classification envisages three main morphologies defined following the channel-reach morphology classifications proposed by Grant et al. (1990) and Montgomery and Buffington (1997): cascade, step pool and rapid with steps. Cascade reaches are characterized by longitudinally and laterally disorganized bed material typically consisting of cobbles and boulders, whereas step-pool reaches are characterized by longitudinal steps formed by large clasts organized into discrete channelspanning accumulations that separate pools containing finer material (Montgomery and Buffington, 1997). The cascade reach at the outlet of the rocky gorge, which is characterised by abundant large boulders, was classified as a boulder cascade to better emphasize the coarse grain size of this channel reach. One of the two channel reaches previously classified as riffle-pool was interpreted as a rapid with steps (Grant et al., 1990), characterised by an organisation of small and medium boulders into irregular ribs oriented perpendicular to the channel and exposed to low flow (Fig. 3), whereas the second displays sufficiently defined steps and was ascribed to two neighbouring step-pool and stepped cascade reaches. The term “mixed reach,” which identifies channel reaches similar to step-pool sequences but with coarse particle bars deposited upstream of the boulder steps or adjacent to isolated big boulders (Billi et al., 1998), was abandoned, and the channel reaches previously classified as “mixed” were www.hydrol-earth-syst-sci.net/14/393/2010/ classified as cascades. The presence of large woody debris in the studied channel is very limited and does not influence channel morphology. The upper and lower ends of channel reaches were checked using GPS, leading to slight differences from the previous sequences of step-pool and mixed reaches. Figure 2 shows a longitudinal profile of the studied channel, with the classification of channel bed morphologies. 3 LiDAR data LiDAR and photographic data were acquired from a helicopter using an ALTM 3100 OPTECH and a Rollei H20 digital camera. The survey was carried out in snow-free conditions in October 2006, flying at an average altitude of 1000 m above ground level. The flying speed was 80 knots, the scan angle 20 degrees and the pulse rate 71 KHz. The survey design point density was specified to be greater than 5 points/m2 , recording up to four returns, including the first and last. The shallow depth of water in the analysed stream at the time of the LiDAR flight (not exceeding a few tens of centimetres) made it possible to survey the channel bed without using a bathymetric LiDAR sensor. However, some attenuation of the LiDAR signal has been observed in the pools of step-pool reaches. LiDAR point measurements were filtered into returns from vegetation and bare ground. A bare ground dataset was generated by eliminating non-ground points by manual editing, through visual inspection of aerial photos and using 2-D plots (elevation versus distance along the channel, Fig. 4), which allows points located in vegetated areas to be identified. This approach, which is fairly time-consuming but still affordable thanks to the small extent of the analysed area (8144 m2 ), was chosen because it permits a careful control on the filtering process. Manual filtering removes Hydrol. Earth Syst. Sci., 14, 393–405, 2010 396 S. Trevisani et al.: Reading the bed morphology of a mountain stream Merwade et al., 2005; Legletier and Kyriakidis, 2007a). The coordinate transformation was performed using spatial tools available in most GIS packages, such as spatial queries, spatial joins and shape discretisations. In the present study, we used ArcGIS 9.3 (ESRI, 2008). The procedure is composed of the following steps: – Derivation of the centreline: this is computed as the line equidistant from the two lines defining the boundary of the domain. This line must be smoothed to avoid having a radius of curvature less than the half-width of the channel; otherwise it would not be possible to perform the transformation (Goff et al., 2004; Legletier and Kyriakidis, 2007b). Fig. 3. Picture of the rapid with steps reach. misclassification errors (i.e., points falling on trees classified as ground), which would have appeared as strong local outliers from a geostatistical perspective. The filtering of less obvious errors, like those related to instrumentation, is handled during the interpolation sequence. The average density of ground points was 5.46 points/m2 , with a standard deviation of 3.69. This value indicates that the manual filtering approach can preserve the high density of the raw dataset. Mean and median separations between ground points, determined with the nearest neighbour technique, are 0.19 and 0.17 m, respectively. 4 4.1 Methods Domain delineation and coordinate transformation This study required the definition of the spatial analysis domain, which corresponds to the channel bed, by delineating the lower edges of the banks. Elevation data falling on the channel banks, which could bias the interpolation and computation of morphological indices, were excluded from the analysis. In order to digitise and georeference in a GIS the outline of the channel bed, different thematic layers were analysed interactively with the support of aerial orthophotos with a horizontal resolution of 0.15 m. The edges corresponding to the feet of the banks were detected using maps of shaded relief, slope, roughness (Cavalli et al., 2008) and openness index (Yokoyama et al., 2002). These maps were derived from a preliminary DTM with pixels of 0.5 m×0.5 m. The derived channel bed outline ranges in width from 1.5 to 10.5 m, with a mean value of 5 m (first and third quartiles of 4 m and 6 m). Subsequent processing of LiDAR data was performed on a transformed coordinate system (Fig. 5). Spatial analysis in transformed coordinate systems has found many applications in the field of geostatistical geomodelling, in environmental geostatistics, and in river geomorphology (Deutsch, 2002; Hydrol. Earth Syst. Sci., 14, 393–405, 2010 – Centreline discretisation: the smoothed centreline is discretised as points with a close spacing (step of 1 cm), ensuring a high-accuracy transformation of the coordinate system. – Computation of stream-centred coordinates: in the new system, the position of a given point in the river bed is defined by an abscissa, RX, and an ordinate, RY. RX indicates the distance (positive or negative depending on the side) from the nearest point discretising the centreline. RY corresponds to the position of the nearest point of discretisation on the centreline (Fig. 5). The derivation of the centreline is quite straightforward in the Rio Cordon, which has a narrow channel with low sinuosity. In the case of high sinuosity and wider streams (e.g., meandering rivers), the centreline derivation could require more trials in the smoothing process, to get a curvature radius less than half of the channel width. The transformation of the coordinate system from a Cartesian to a stream-centred one has advantages both for the interpolation (i.e., DTM generation) and for spatial analysis. DTM interpolation on the stream-centred coordinate system allows the various parameters and models (variogram anisotropies, search ellipse, trend calculation) to be adjusted in a system in which the main directions, along and transverse to the stream flow direction, have a sound physical basis. In the transformed coordinate system, the calculation of directional morphological indices in directions orthogonal and parallel to river flow becomes straightforward, and indices evaluating the asymmetry of bedform organisation with respect to the centreline can be easily calculated. Moreover, the stream-centred coordinate system emulates the viewpoint of the expert in the field, who surveys a river channel taking as the main directions the flow direction of water (channel profile) and the direction orthogonal to it (cross-sections). 4.2 DTM derivation The DTM was derived via ordinary block kriging interpolation, with a pixel size of 0.5 m. The DTM resolution was chosen taking into account the spatial density of LiDAR data www.hydrol-earth-syst-sci.net/14/393/2010/ S. Trevisani et al.: Reading the bed morphology of a mountain stream 397 Fig. 4. Elevation of LiDAR points in a sample area of the study domain, plotted along the centreline of the channel (a). Non-ground points (highlighted in a and in b), which correspond to trees, are easily recognizable, and can be expunged from the dataset. and the need to describe local (i.e., large-scale) morphological features. The high data density (5.46 points/m2 ) allowed us to use a simplified interpolation approach, avoiding the use of more complex procedures, such as universal kriging or local kriging (Stroet et al., 2005). Block kriging interpolation was carried out using an elliptical search window with a semiaxis of 1 m in the flow direction and of 2 m across it. The short semiaxis in the direction of flow was dictated by the need to attenuate the effect of trend, i.e., the downstream decrease of elevation. For the distances of interest, the variogram was modelled by a spherical variogram with a range of 4 m and a sill of 0.3 m2 plus a nugget of 0.01 m2 . We considered the nugget as representative of measurement error and not of microvariability (Cressie, 1993), corresponding to white noise with a standard deviation 0.1 m. The imposition of the nugget as error permits us to filter out the white noise (i.e., not spatially correlated) component of error measurement. Then, the interpolation via block kriging, filtering out the intra-block spatial variability, is capable of removing short-wavelength spatially-correlated error. The quality of the derived DTM was evaluated quantitatively by means of standard deviation maps, cross-validation results and maps, and qualitatively during the field surveys visually comparing the shaded relief map with the real streambed morphology. 4.3 Morphometric indices Morphometric indices have been computed to characterise the morphology of the channel bed of the Rio Cordon. Except for the local channel slope, which is calculated based on elevation data, the other indices were calculated based on the residuals of elevation (Wooldridge and Hickin, 2002; Zimmermann et al., 2008; Trevisani et al., 2009), which were derived by removing the large-scale spatial variation from www.hydrol-earth-syst-sci.net/14/393/2010/ the interpolated DTM, i.e., removing the trend. The surface of residuals represents the high-frequency spatial variability component of the streambed morphology. The trend surface was derived by means of a local polynomial approach, fitting a planar model in a moving window covering the whole width of the stream for a length of 10 m in the direction of flow. The morphological indices were calculated within moving windows; the calculated value was attributed to the centre of the window (Fig. 6). The indices were computed every 1 m along the stream centreline, using moving windows covering the whole width of the channel. Both the DTM resolution and the moving window sizes were chosen taking into account the sizes of the topographic features being studied. The DTM cell size of 0.5 m is suitable to correctly represent the step pools’ grain size, which is characterised in the Rio Cordon by a D84 of about 1 m (Lenzi, 2001). Two different moving windows lengths (2 and 10 m) along the flow direction were used. The former (SW, small window) was chosen to detect single morphological units or small-scale features like steps, small pools and boulders; the latter (LW, large window) can describe the average spatial characteristics of a longer stream reach covering more morphological units (the mean step-pool spacing reported by Lenzi (2001) ranges from about 2.5 to 5 m), thus facilitating the recognition of morphological sequences. The computation of the indices at regular intervals along the longitudinal channel profile, for windows covering the whole channel width, makes possible to revert the 2-D analysis to a 1-D issue, which can be tackled using time series analysis (Malamud and Turcotte, 1999); the distance along the stream centreline corresponding to time. An important advantage arises in the interpretation of the indices, which can be plotted as a profile along the channel, making it straightforward to visually compare with the topographic channel profile and morphological classification. Hydrol. Earth Syst. Sci., 14, 393–405, 2010 398 S. Trevisani et al.: Reading the bed morphology of a mountain stream (a) (b) Fig. 6. (a) the SW (blue) and LW (red) moving windows are superposed on the surface of residuals. (b) a sketch of the computation of directional variograms via SW windows is presented. Fig. 5. Coordinate transformation from a geographical (a) to a stream-centred system (b). A number of trials on different indices and computation parameters have been carried out to select the indicators to analyse the morphology of the Rio Cordon. The tested indices span from terrain analysis indices, such as directional slope and directional curvature, to spatial statistical indices, such as directional variograms calculated at different lags and parameters of variogram models fitted on experimental local variograms (Butler et al., 2001). The selection of geomorphometric indices has been based on a balance of complexity and informative content. One basic morphometric index is the longitudinal slope, calculated by fitting a line, via least squares, to the elevation versus RY-coordinate plot. The regression was carried out on the subset of points falling within each moving window. This index, when calculated within SW (SW slope), can detect local features of channel bed morphology. When calculated on LW (LW slope), it gives overall information on the hydraulic energy of the investigated reaches. It is important to stress that the slope index is not calculated along the centreline of the channel, but considers the whole channel width, and is thus able to synthesise the average slope of the 2-D channel on a 1-D representation. The effectiveness of a 1-D representation considering the whole channel width for morphological analysis on step pool channels was demonstrated in a flume study by Zimmermann et al. (2008). The authors used average long profiles, constructed by averaging the bed elevation across the channel for different cross-sections, to objectively classify step pool units. Hydrol. Earth Syst. Sci., 14, 393–405, 2010 A step toward a finer characterisation of channel-bed morphology is the calculation of directional variograms γ (h) (1) for the shortest lag (0.5 m) (Fig. 6). γ (h) = X 1 N(h) [z(uα ) − z(uα +h)]2 2N (h) α=1 (1) with γ (h) = γ (−h) As outlined in Eq. (1), the variogram gives, for each value of the vector h, an estimate of the half-mean squared differences between the N couples of points z(uα ) separated by the vector h. Thus, the variogram, for a specific lag, gives a measure of spatial variability functional to separation distance, i.e., a scale-dependent directional roughness index. The variogram was calculated along the RX and RY directions. The variogram calculated along RY (Fig. 6) can be viewed as a proxy for the hydraulic resistance (higher variogram values indicate higher flow resistance). The variogram calculated in the directions orthogonal to the river flow (RX) is intended to give information on the transverse organisation of the sediment. The calculation of the indices was performed via ad hoc routines developed inside the open-source statistical programming and graphics environment R (R Development Core Team, 2009). For variogram calculation, some of the functions of the Gstat R packages were used (Pebesma et al., 2004). A great improvement in the readability of the directional variograms index can be achieved by calculating local anomalies (A) as the difference between the variograms www.hydrol-earth-syst-sci.net/14/393/2010/ S. Trevisani et al.: Reading the bed morphology of a mountain stream calculated in small and large windows (a method analogous to the sharpening technique applied in image processing): A = γSW − γLW (2) In Eq. (2), the value of the variogram on LW is calculated from the original data, but similar results can be obtained (the differences are related to the different cardinalities of data encountered inside the moving windows) by calculating it, by means of moving averages, directly from the time series of the variogram calculated on SW, which represents the signal at the highest possible resolution. The interpretation of the local anomaly is straightforward: it determines whether the local variability on SW is higher or lower than the mean variability calculated on LW. This approach can be generalised to a multiresolution analysis by increasing the size of the windows. If S(tn ) is the discrete time series representing the analysed signal (here the variogram calculated on SW for a lag of 0.5 m), defined at times tn (tn = nδ; n=1, 2, 3. . . N; where δ is the time step, and N the number of time steps), the local anomaly at the level of resolution l, Al (tn ), can be derived as: Al (tn ) = m(λl ,S(tn )) − m(λl+1 ,S(tn )), (3) where m(λl ,S(tn )) indicates the moving average process of the signal S(tn ) centred at every time step tn within moving windows of size λl , where λl = λ1 ,λ2 ,λ3 ,. . . , λL indicate windows of increasing size up to the level of resolution L. According to this formulation, the anomalies are calculated at different scales: by increasing the size of the windows it is possible to analyse the characteristics of the signal at lower frequencies. This approach gives similar results, at least for the datasets analysed in our case studies, to an analysis performed via wavelets (Lark and Webster, 1999; Percival and Walden, 2000), and, thanks to its simplicity, could be preferred to wavelets for an explorative analysis of data. For other tasks, such as image segmentation, signal compression and pattern recognition, wavelets analysis, after a careful choice of wavelet type and calculation parameters, is necessary. In this case, to get a more smoothed signal, we used a modified version of Eq. (3), practically calculating the means at a given level l from the means of level l − 1 and not from the original signal: Al (tn ) = m(λl ,ml−1 (tn )) − m(λl+1 ,ml (tn )) (4) with m0 (tn ) = S(tn ) 4.4 Field surveys Field surveys were carried out with two main objectives: – to revise the previous classification of channel bed morphology (Billi et al., 1998, Lenzi, 2001, Cavalli et al., 2008) for achieving a more detailed identification of morphological sequences; www.hydrol-earth-syst-sci.net/14/393/2010/ 399 – to check the correspondence between the results of digital terrain analysis and ground truth. The revised classification of channel-bed morphology of the Rio Cordon has been presented in the description of the study area (Figs. 2 and 3). The field survey aimed at evaluating the capability of the morphometric indices to depict morphological features of the channel bed was conducted using mobile GIS software (ESRI ArcPad 7.1) installed on a rugged tablet PC with GPS. The thematic layers of the GIS were as follows: – Orthophoto at a high resolution (0.15 m); – Shaded relief maps of the DTM; – Residuals of the interpolated surface; – Contour lines. Longitudinal profiles of the calculated indices (directional variograms and local anomalies) were used for a field comparison with morphological features of the channel. The mobile GIS integrated with GPS made it possible to identify in the field the exact position along the centreline in the transformed coordinate system. The survey of the main channel of the Rio Cordon allowed us to compare values of the proposed morphometric indices with actual channelbed morphology. Furthermore, the representativeness of the LiDAR-derived DTM was visually checked during field surveys, paying particular attention to the areas where a deterioration of LiDAR data quality or density could be expected (i.e., vegetation and pools with relatively deep water). 5 Results Figure 7 presents the slope calculated on SW (Fig. 7a) and LW (Fig. 7b) along the channel profile. The patterns on the graphs of the SW and LW slopes are similar, with the strongest oscillations of slopes located in the boulder cascade reach and the lowest at the rapid with step reach. Obviously, the spatial variability of slope and its range of variations are larger for SW, due to its capability to detect high detail features. Negative values are visible in the plot of SW slope. Negative slopes in step pool reaches correspond to a rise of the channel bed at the downstream end of pools; in the boulder cascade reach and the cascade reaches, these are mainly related to the presence of scattered large boulders. The slope index is useful for a general characterisation of the channel profile and should be used preliminarily to the geostatistical indices, which focus more closely on channel bed roughness. The directional variograms calculated along the RX and RY directions for a lag of 0.5 m are presented for SW (SWG) and LW (LWG) in Fig. 8. Hydrol. Earth Syst. Sci., 14, 393–405, 2010 400 S. Trevisani et al.: Reading the bed morphology of a mountain stream Fig. 7. Longitudinal slope [m/m] calculated on SW and LW, plotted along the centreline of the streambed. The overall patterns of the RX and RY variogram graphs are similar, although the RX variogram is generally more scattered, especially in the lower part of the channel. This could be related to the fact that the transverse organisation of the sediments, also given the limited width of the stream, is not clearly distinctive for the different morphologies. Accordingly, the discussion will largely focus on the analysis of variogram-based indices calculated along RY (Fig. 8c and d). – In the upper part of the profile (from 0 to 250 m), which mostly corresponds to the boulder cascade reach, the variability is very high. The highest values are observed in the upper 150 m, where the Rio Cordon flows between steep rocky slopes. Here large boulders accumulated by rockfalls cause a high roughness of the channelbed surface. – From 250 to 480 m, the roughness of the step pool and cascade reaches is markedly lower than in the upstream boulder cascade reach. – In correspondence to the rapid with steps reach and part of the cascades reaches from 480 to 630 m, the variability is very low. – From 630 to 750 m, the variability increases again, with the presence of definite spikes (cascade and step-pool reaches). – From 750 to 980 m, the graph appears less variable (in the same zone the RX variogram graph is still quite erratic); this happens in correspondence to a quite long cascade reach. – In the remaining part of the stream profile, a succession of high and low variability zones is present, with Hydrol. Earth Syst. Sci., 14, 393–405, 2010 clear spikes related to the step-pool sequences, and to the presence of large boulders. In this section, an alternating sequence of step-pool and cascade reaches is present, with a morphology strongly influenced by lateral tributaries and the deposits of debris-flows. Figure 9 presents multiscale local anomalies (Eq. 4) of RY together with the original signal (bottom, Fig. 9a), represented by the SW variogram along RY (SWGry ), and the remaining smoothed signal (Fig. 9f), in a format analogous to the one used in wavelet analysis. The main patterns visible in the anomaly graphs persist along the first three levels of resolution, while in the fourth level the structure is quite simplified. The smoothed signal (Fig. 9f, moving window of 80 m), which can be interpreted as a simplification of the original signal (Fig. 9a), outlines the main patterns just described for SWGry . In the present study because of our interest in characterising the morphology at the finest possible scale, we are especially interested in the anomalies at the highest level of resolution (Fig. 9b, anomaly A1). Even if the interpretation of the signal is similar to that of the variogram on RY computed on SW in the anomaly plot, the zones of different variability are more distinct. Peculiar patterns and values observed in the plot of anomalies at high resolution, in combination with slope and variogram indices, have been checked in the field, looking for corresponding features in channel bed morphology (Fig. 10). Figure 10a depicts a detail of the upstream part of Rio Cordon, classified as a boulder cascade, which is characterised by a scattered pattern with a wide range of variation and the presence of marked positive spikes. In Fig. 10b, two moderate negative peaks separated by a strong positive peak represent, respectively, the smoothed shapes of two pools and the high variability of the step. Figure 10c shows a pattern of anomalies characterised by small variations and by a quite www.hydrol-earth-syst-sci.net/14/393/2010/ S. Trevisani et al.: Reading the bed morphology of a mountain stream 401 Fig. 8. RX and RY directional variogram values calculated on SW and LW for lag 0.5 m. symmetric distribution of positive and negative anomalies in the reach classified as a rapid with steps. Here steps built by relatively small grain size material are very low and short: this explains the plane trend in the anomaly graph. The highest positive peaks in the lower part of the channel are indeed related to significant structural roughness features, such as boulder clusters and big steps (Fig. 10d and e). In Fig. 10d, the isolated strong positive spike represents the big step visible in the photo. In Fig. 10e, an isolated step pool appears as negative peaks divided by a strong positive peak. The pattern and range of variation of the anomalies represent the different character of the analysed stream portion and can detect particular morphological units, like step pool sequences shown as negative peaks divided by a strong positive peak. In the Rio Cordon, the pools, highlighted by negative anomalies, can be considered as the more distinctive features between the different channel bed morphologies, whereas clear positive anomalies can be found in all morphologies, but rapid with steps. The field observations summarised above show that the index is also quite robust in the presence of an articulate morphology. A previous study on the channel bed of the Rio Cordon (Cavalli et al., 2008) analysed the main morphological www.hydrol-earth-syst-sci.net/14/393/2010/ classes by means of morphological indices of roughness based on LiDAR-derived DTM. The analysis was conducted by pooling all channel reaches belonging to the same morphological class. The comparison between different morphological classes demonstrated that roughness indices of the channel bed surface have significantly different values for different morphologies. In the present study, a geostatistical approach has permitted an enhanced representation of the spatial variability of channel bed morphology, allowing the recognition of even single morphological units. 6 Conclusions The main outcomes of the geomorphometric analysis of channel bed morphology of the Rio Cordon are reported below. – Working in a river-centred coordinate system was very useful for both the interpolation procedure and the geomorphometric study. The transformed coordinate system makes is feasible to calculate directional geomorphometric indices in the main process-based directions and permits us to synthesise the 2-D spatial information Hydrol. Earth Syst. Sci., 14, 393–405, 2010 402 S. Trevisani et al.: Reading the bed morphology of a mountain stream Fig. 9. Anomalies of the RY variogram for a lag of 0.5 m, for four levels of spatial resolution. The results are presented in analogy to wavelets analysis. (a) Original signal (i.e., longitudinal variogram for a lag of 0.5 m). (b–e) Multiscale anomalies for decreasing levels of resolution calculated according to Eq. (4). The sizes of the windows are indicated in the figure. (f) Remaining smoothed signal. according to a 1-D data representation, easily interpretable by means of various graphical, statistical and analytical tools. – A slope index computed for the whole width of the channel bed using moving windows of different sizes makes it possible to depict the variability of slope along the channel profile, including the presence of negative slopes, mostly associated with pools in step-pool channel reaches. Hydrol. Earth Syst. Sci., 14, 393–405, 2010 – The variogram computed in the direction of flow on the residuals of the DTM clearly shows changes in topographic variability along the streambed profile that depend on the different morphologies of the channel bed. Plots of local anomalies, computed as the difference between directional variograms calculated in small and large moving windows, enhance the highdetail variations in the topographic roughness of the channel bed. The field check has shown a satisfactory www.hydrol-earth-syst-sci.net/14/393/2010/ S. Trevisani et al.: Reading the bed morphology of a mountain stream 403 Fig. 10. Comparison of local anomalies of the RY variogram for a lag of 0.5 m with field observations of the channel bed morphology for some portions (highlighted by boxes along the profile). Boulder cascade (a), step-pool sequence (b); rapid with step (c); isolated big step (d); and isolated step pool (e). www.hydrol-earth-syst-sci.net/14/393/2010/ Hydrol. Earth Syst. Sci., 14, 393–405, 2010 404 S. Trevisani et al.: Reading the bed morphology of a mountain stream correspondence between values of local anomalies and particular features of channel morphologies, such as pools, large steps, clusters of boulders and pools. These morphological units, although also visible in topographic profiles from high-resolution data (Cavalli et al., 2008), need directional indices of topographic variability, such as directional variograms and local anomalies, to be adequately emphasised. – Boulder cascades and rapids with step reaches can be easily recognised in graphs of longitudinal variograms and local anomalies, whereas other classes, such as step pools and cascades, show similar general patterns. Nevertheless, well-defined pools between steps which represent distinctive features between step pool and cascade morphologies, can be easily recognised as negative values in the graph of local anomalies (Fig. 10b and e). The indices proposed in this paper are thus effective for recognising single morphological units (Fig. 10), whereas further procedures are needed for the automatic or semi-automatic classification of the morphology at the channel reach scale. The Rio Cordon is a challenging case for the geomorphometric study of streambeds because of the presence of many short channel reaches with different characteristics (Fig. 2) and fuzzy transitions between the various morphologies. Moreover, the narrow channel width and DTM resolution limit the calculation of variogram to the shortest lag, leading to results similar to those achievable by means of directional derivatives. As a consequence, the potentials of the proposed geostatistical indices are not fully exploited. However, directional derivatives are less general and scalable than the variogram (or other spatial continuity indices such as the more robust madogram), which could be fully exploited in other morphological contexts (e.g., larger streambeds) and with higher resolution datasets acquired, for example, using terrestrial LiDAR. In many alpine headwater streams the morphological conditions of the channel are complex, like those of the Rio Cordon, and topographic datasets from aerial LiDAR, similar to that used in this study, are increasingly available. For these reasons, the Rio Cordon can be deemed as representative of achievements and problems that may occur in the geomorphometric analysis of channel bed morphology of small alpine streams by means of LiDAR data. Acknowledgements. The authors wish to thank Francesco Comiti and Alessandro Vianello for their advice in the revision of the classification of channel bed morphology of the Rio Cordon and the referees for their useful comment. Edited by: P. Molnar Hydrol. Earth Syst. Sci., 14, 393–405, 2010 References Ardizzone, F., Cardinali, M., Galli, M., Guzzetti, F., and Reichenbach, P.: Identification and mapping of recent rainfall-induced landslides using elevation data collected by airborne Lidar, Nat. Hazards Earth Syst. Sci., 7, 637–650, 2007, http://www.nat-hazards-earth-syst-sci.net/7/637/2007/. Billi, P., D’Agostino, V., Lenzi, M. A., and Marchi, L.: Bedload, slope and channel processes in a high-altitude alpine torrent, in: Gravel-bed Rivers in the Environment, edited by: Klingeman, P. C., Beschta, R. L., Komar, P. D., and Bradley, J. B., Water Resources Publication, LLC, Highland Ranch, Colorado, USA, 15–38, 1998. Brown, A. G., Carey, C., Erkens, G., Fuchs, M., Hoffmann, T., Macaire, J., Moldenhauer, K., and Walling, D.: From sedimentary records to sediment budgets: Multiple approaches to catchment sediment flux, Geomorphology, 108(1–2), 35–47, 2009. Butler, J. B., Lane, S. N., and Chandler, J. H.: Characterization of the Structure of Stream-Bed Gravels Using Two-Dimensional Fractal Analysis, Math. Geol., 33(3), 301–329, 2001. Cavalli, M., Tarolli, P., Marchi, L., and Dalla Fontana, G.: The effectiveness of airborne LiDAR data in the recognition of channelbed morphology, Catena, 73(3), 249–260, 2008. Charlton, M. E., Large, A. R. G., and Fuller, I. C.: Application of airborne LiDAR in river environments: the River Coquet, Northumberland, UK, Earth Surf. Proc. Land., 28(3), 299–306, 2003. Chin, A.: The morphological structure of step–pools in mountain streams, Geomorphology, 27(3–4), 191–204, 1999. Cressie, N. E.: Statistic for Spatial Data, revised edition, John Wiley & Sons Inc., New York, 900 pp., 1993. Deutsch, C. V.: Geostatistical Reservoir Modeling, Oxford University Press, New York, 376 pp., 2002. ESRI: ArcGIS version 9.3. Redlands, California, USA, 2008. Glenn, N. F., Streutker, D. R., Chadwick, D. J., Thackray, G. D., and Dorsch, S. J.: Analysis of LiDAR-derived topographic information for characterizing and differentiating landslide morphology and activity, Geomorphology , 73(1–2), 131–148, 2006. Goff, J. A. and Nordfjord, S.: Interpolation of fluvial morphology using channel-oriented coordinate transformation: a case study from the New Jersey shelf, Math. Geol., 36(6), 643–658, 2004. Goovaerts, P.: Geostatistics for Natural Resources Evaluation, Oxford University Press, New York, 481 pp., 1997. Grant, G. E., Swanson, F. J., and Wolman, M. G.: Pattern and origin of stepped-bed morphology in high-gradient streams, Western Cascades, Geol. Soc. Am. Bull., 102(3), 340–352, 1990. Hengl, T. and Reuter, H. I.: Geomorphometry: Concepts, Software, Applications. Developments in Soil Science, Elsevier, Amsterdam, 33, 772 pp., 2009. Höfle, B., Vetter, M., Pfeifer, N., Mandlburger, G., and Stötter, J.: Water surface mapping from airborne laser scanning using signal intensity and elevation data, Earth Surf. Proc. Land., 34(12), 1635–1649, 2009. Jones, A. F., Brewer, P. A., Johnstone, E., and Macklin, M. G.: High-resolution interpretative geomorphological mapping of river valley environments using airborne LiDAR data, Earth Surf. Proc. Land., 32(10), 1574–1592, 2007. Lark, R. M. and Webster, R.: Analysis and elucidation of soil variation using wavelets, Eur. J. Soil Sci., 50, 185–206, 1999. www.hydrol-earth-syst-sci.net/14/393/2010/ S. Trevisani et al.: Reading the bed morphology of a mountain stream Legletier, C. J. and Kyriakidis, C.: Spatial prediction of river channel topography by kriging, Earth Surf. Proc. Land., 33(6), 841– 867, 2007a. Legletier, C. J. and Kyriakidis, C.: Forward and Inverse Transformation between Cartesian and Channel-fitted Coordinate Systems for Meandering Rivers, Math. Geol., 38(8), 927–957, 2007b. Lenzi, M. A.: Step-pool evolution in the Rio Cordon, Northeastern Italy, Earth Surf. Proc. Land., 26(9), 991–1008, 2001. Magirl, C. S., Webb, R. H., and Griffiths, P. G.: Changes in the water surface profile of the Colorado River in Grand Canyon, Arizona, between 1923 and 2000, Water Resour. Res., 41, W05021, doi:10.1029/2003WR002519, 2005. Malamud, B. D. and Turcotte, D. L.: Self affine time series: 1. Generation and analyses, Adv. Geophys., 40, 1–90, 1999. Mason, D. C., Scott, T. R., and Wang, H.-J.: Extraction of tidal channel networks from airborne scanning laser altimetry, ISPRS J. Photogramm., 61(2), 67–83, 2006. Merwade, M. V., Maidment, D. R., and Hodges, B. R.: Geospatial Representation of River Channels, J. Hydrol. Eng., 3, 243–251, 2005. Milzow, C., Molnar, P., McArdell, B. W., and Burlando, P.: Spatial organization in the step-pool structure of a steep mountain stream (Vogelbach, Switzerland), Water Resour. Res., 42, W04418, doi:10.1029/2004WR003870, 2006. Montgomery, D. and Buffington, J. M.: Channel-reach morphology in mountain drainage basins, Geol. Soc. Am. Bulletin, 109(5), 596–611, 1997. Pebesma, E. J.: Multivariable geostatistics in S: the gstat package, Comput. Geosci., 30(7), 683–691, 2004. Percival, D. B. and Walden, A. T.: Wavelets Methods for Time Series Analysis, Cambridge University Press, New York, 620 pp., 2000. Pike, R. J.: Geomorphometry – diversity in quantitative surface analysis, Prog. Phys. Geog., 24(1), 1–20, 2000. R Development Core Team R: A language and environment for statistical computing, R Foundation for Statistical Computing, Vienna, Austria, 2009. www.hydrol-earth-syst-sci.net/14/393/2010/ 405 Scheidl, C., Rickenmann, D., and Chiari, M.: The use of airborne LiDAR data for the analysis of debris flow events in Switzerland, Nat. Hazards Earth Syst. Sci., 8, 1113–1127, 2008, http://www.nat-hazards-earth-syst-sci.net/8/1113/2008/. Stroet, C. B. M. and Snepvangers, J. J. J. C.: Mapping Curvilinear Structures with Local Anisotropy Kriging, Math. Geol., 37(6), 635–649, 2005. Tarolli, P. and Dalla Fontana, G.: Hillslope-to-valley transition morphology: new opportunities from high resolution DTMs, Geomorphology, 113(1–2), 47–56, 2009. Thoma, D. P., Gupta, S. C., Bauer, M. E., and Kirchoff, C. E.: Airborne laser scanning for riverbank erosion assessment, Remote Sens. Environ., 95(4), 493–501, 2005. Trevisani, S., Cavalli, M., and Marchi, L.: Variogram maps from LiDAR data as fingerprints of surface morphology on scree slopes, Nat. Hazards Earth Syst. Sci., 9, 129–133, 2009, http://www.nat-hazards-earth-syst-sci.net/9/129/2009/. Vianello, A., Cavalli, M., and Tarolli, P.: LiDAR-derived slopes for headwater channel network analysis, Catena, 76(2), 97–106, 2009. Wilson, J. P. and Gallant, J. C.: Terrain Analysis: Principles and Applications, John Wiley & Sons INC., New York, 479 pp., 2000. Wohl, E. E., Madsen S., and MacDonald, L.: Characteristics of log and clast bed-steps in step–pool streams of northwestern Montana, USA, Geomorphology 20(1–2), 1–10, 1997. Wooldridge C. L. and Hickin, E. J.: Step-pool and cascade morphology, Mosquito Creek, British Columbia: a test of four analytical techniques, Can. J. Earth Sci., 39(4), 493–503. 2002. Yokoyama, R., Shirasawa, M., and Pike, R.: Visualizing topography by openness: a new application of image processing to digital elevation models, Photogramm. Eng. Rem. S., 68(3), 257–266. 2002. Zimmermann, A. E., Church, M., and Hassan, M. A.: Identification of steps and pools from stream longitudinal profile data, Geomorphology, 102(3–4), 395–406, 2008. Hydrol. Earth Syst. Sci., 14, 393–405, 2010