Skip to navigation – Site map

HomeNumérosvol. 11 - n° 4Measurement of DEM roughness usin...

Measurement of DEM roughness using the local fractal dimension

Mesure de la rugosité des MNT à l’aide de la dimension fractale
Hind Taud and Jean-François Parrot
p. 327-338

Abstracts

The relationships between geological features and DEM surface roughness are studied using fractal geometry analysis. The estimate of the fractal dimension in a 3D space is performed locally on DEM surfaces using an adaptive “box-counting” technique. This procedure has been applied to two regions chosen to represent differences in lithological and tectonic conditions. The first one corresponds to a faulted, homoclinal sedimentary sequence and the second to a compound stratovolcano. In both cases, the results obtained show that the local fractal dimension detects distinct morphologic features. In the first case (Vittel, eastern France), the local fractal dimension detects the limits of the geological units as borders of homogeneous zones or transitions between different dissection depths. Fine structural lines underlying the Vittel fault can also be extracted. In the case of Mount Ararat (Eastern Turkey), different classes can be distinguished within the volcanic formations using a statistical analysis of the fractal dimension values. The treatment can also provide information about the local and regional geological structures. These first results show that measuring DEM surface roughness represents a helpful tool for extracting and mapping morphometric features.

Top of page

Errata

Article reçu le 20 septembre 2004, accepté le 24 août 2005

Editor's notes

Acknowledgements
We thank the reviewers for their helpful criticism and suggestions, the Mexican Petroleum Institute (IMP) and the Geographical Institute (UNAM) for their support.

Full text

Introduction

1As a Digital Elevation Model (DEM) is a representation of a surface, several algorithms have been developed to study the surface properties of DEMs and to provide a large set of descriptor attributes (Wilson and Gallant, 2000). Given that altitudes correspond to grey levels, it is then possible to describe, quantify, and model rough surfaces using image-processing techniques in the spatial and frequency domains and to apply pattern recognition processes. Textural analysis is closely linked to roughness assessment. The concept of texture is quite difficult to define as it often includes the notions of roughness, regularity, contrast, and thinness. Numerous authors have attempted to outline and clarify this notion. J. C. Russ (1999) proposed that roughness equates with the high frequencies of the 2D signal constituting a DEM. He pointed out that a DEM contains three levels of information. The first level, related to the lower frequencies, coincides with its form; the second level, related to the middle frequencies, corresponds to its waviness; and the third level corresponds to its roughness.

2Some recent techniques in image analysis use fractal and/or multi-fractal approaches to characterize the texture of a greyscale image or the roughness of a surface. Concerning DEM surfaces, several authors have demonstrated that the topography of the Earth generally exhibits fractal characteristics and that the relief preserves the same statistical characteristics over a wide range of scales (Huang and Turcotte, 1989; Klinkenberg and Goodchild, 1992). By describing the terrain as a fractal surface, the local or the global analysis of DEM roughness reveals interpolation artefacts, provides an assessment of the DEM quality (Polidori et al., 1991; Datcu et al., 1996), and can assist in understanding erosional phenomena (Chase, 1992; Cheng et al., 1999).

3As every landscape appears to have a particular fractal dimension value, it implies that this dimension calculated locally has to be related with the particular features of the local landscape units. On the other hand, each local landscape unit presents its own morphologic feature in relation to the erosion processes and mainly to the nature of the basement. Thus, the purpose of this research is to study the relationships between geomorphic features and the surface roughness of a DEM by locally measuring fractal dimension in 3D space.

4The purpose here does not include demonstrating the fractal behaviour of the DEM nor comparing the performances of different fractal surface estimators. The first section of the paper presents a brief overview of fractal geometry and explains the notion of local fractal dimension. The second section describes the procedure, while the third section reports the results obtained in two case studies. These training areas have been chosen according to diverse lithological conditions. The well known geological features of the sedimentary Vittel region leads to investigate the ability of the treatment to extract these features in accordance to their morphology. Taking into account the former results, the treatments applied to the Mount Ararat illustrate the reliability of the method to emphasize local roughness differences related to regional tectonism and to detect peculiar and unknown features that can be useful for geological and geomorphological mapping.

Fractal geometry and local fractal dimension

Fractal geometry

5Fractal geometry, introduced and developed by B. Mandelbrot (1982), provides a mathematical description of a wide range of natural forms and phenomena. Fractal objects are defined as scale-invariant (self-similar or self-affine). This means that the fractal object can be presented as an assemblage of rescaled copies of itself. Self-similarity occurs when the rescaling is isotropic or uniform in all directions, and self-affinity occurs when the rescaling is either anisotropic or dependent on direction.

6Fractal objects exhibit details at arbitrarily small scales, and they are too complex to be represented in a Euclidean space. Also known as the Hausdorff-Besicovitch dimension (Falconer, 1990), the fractal dimension differs from the more familiar Cartesian or topological dimension. In this last case, integer values are required: 1 for a line, 2 for a surface, and 3 for a volume. The fractal dimension measures the complexity of the object. A shape with a higher fractal dimension is more complicated or irregular than one with a lower dimension. For example, a shape with a fractal dimension falling within the range [1, 2] fills more space than a one-dimensional curve and less space than a two-dimensional surface. Self-similarity is defined statistically when it cannot be tested through an infinite range of scales. The statistical fractal behaviour is then related to a given scale range. When the fractal dimension of a particular pattern changes within consecutive ranges of scale, one generally refers to the notion of multi-fractality.

7Several methods are available for estimating the fractal dimension of surfaces, such as the fractional Brownian model (Mark and Aronson, 1984), triangular prism areas (Clarke, 1986), box-counting (Falconer, 1990), and the projective covering method (Xie and Wang, 1999). Measurement may either be direct, when it is applied to grey tones, or indirect as in the case of the examination of “profiles” (one-dimensional transects) or isolines (conversion of the surface into a contour map). In this study, the box-counting method was used because it can be applied to various sets of any dimension and patterns with or without self-similarity (Peitgen et al., 1992). According to K. Falconer (1990) who has discussed the mathematical aspect of the box-counting method, the fractal dimension D can be derived from the relation:

81=NssD or D—logNs/log(s) where Ns corresponds to the number of boxes of size (s) needed to cover the structure. The D value is calculated with the following formula:

Local fractal dimension

9The local measurement of fractal dimensions has been developed to investigate image texture and surface roughness. When this type of measurement is applied to an image I, the results are reported in an image R where the value of each point (i,j) corresponds to its fractal dimension FD calculated in a window of size (2w + 1) x (2w + 1):

10Textural features derived from the local measurement of the fractal dimension have been used for segmentation and classification purposes. Among these techniques, the Fourier power spectrum (Pentland, 1984), the blanket method (Dellepiane et al., 1991), differential box-counting (Sarkar and Chaudhuri, 1994; Chaudhuri and Sarkar, 1995), and the fractional Brownian motion (Chen et al., 1989; Toennies and Schnabel, 1994) are the most common. Using a wavelet technique, M. Datcu et al. (1996) have estimated the local roughness of various DEMs and provided some examples showing that the fractal analysis of these DEMs allows different roughness classes to be distinguished and some artefacts, due to the computation of elevation data, to be detected. For computing the local self-similar properties of a DEM, Y.C. Cheng et al. (1999) have developed a 3D box-counting method applying the triangular prism surface method proposed by K.C. Clarke (1986). They observe that fractal dimension values vary as a function of altitude and have interpreted this phenomenon as reflecting spatial variability in erosional potential. This last method differs from the previous ones because the fractal dimension is estimated in unit areas and not at each point of the image. In our approach, the fractal dimension of each pixel of the DEM surface is calculated inside a moving window centred on this pixel, by using a 3D box-counting adaptive method. The advantage of this treatment is that the volume corresponding to DEM section observed locally is directly related to the altitude values of the pixels defining the DEM surface. The treatment does not require any modelling or interpolation of the surface in order to calculate this dimension.

Procedure

11Inside a cube of size s x s x s centred on the pixel of the studied DEM, the volume corresponding to the surface is defined by a set of voxels. A voxel is an elementary cube, the sides of which are equal to the pixel size. Thus, in the (x,y) space, each point (i,j) of the cube contains a number vs(i,j) of voxels falling between 0 and s. The total number of voxels describing the volume is equal to:

12It is easier to calculate vs (i,j) in a bi-dimensional space as follows:


13where I is the original image, Ic the value of the central pixel, ps the pixel size, h a coefficient that defines the vertical resolution, and s the cube size. The cube is partitioned into boxes of size q varying between 1 and s/2 (fig. 1). Each of these boxes is considered as filled if at least one voxel is contained in this box. In other words, the maximum of vs(i,j)/q in each cell is determined. The computation is done as follows:

14where Max is a function calculating the maximum. The fractal dimension corresponds to the inverse of the slope P=ln(q)/ln(Ns), where q is the size of the box and Ns the total number of filled boxes. When calculating the slope, the coefficient of determination, R2, is computed. The estimate of the fractal dimension depends on various factors as discussed below. As the studied central point is by definition located in the middle of the testing cube, an isolated point, surrounded by values lower than the base of the cube s x s x s, generates a vertical line and a line generates a vertical plane. The slopes are respectively expressed by the following relations: y = -x, y = -2x with R2 = 1. In the same way, a horizontal surface fills the half of the testing cube and the corresponding FD = 3 with R2 = 1 (fig. 2). When the volume presents irregularities its fractal dimension decreases.

Fig. 1 – Local fractal calculation by using a 3-dimensional box counting.
Fig. 1 – Calcul de la dimension fractale locale par la technique tridimensionnelle du « comptage de boîtes ».

Fig. 1 – Local fractal calculation by using a 3-dimensional box counting. Fig. 1 – Calcul de la dimension fractale locale par la technique tridimensionnelle du « comptage de boîtes ».

Example with s= 12 and different grid size q. A: initial topography equivalent to q= 1; B: q= 2; C: q= 3; D: q= 6.
Exemple avec s = 12 et différentes tailles de maille q. A: topographie initiale équivalente à q = 1; B: q = 2; C: q = 3; D: q = 6.

Fig. 2 – Local 3D fractal dimension of a line and a plan with s = 12 and q = 1, 2, 3, 6.
Fig. 2 – Dimension fractale locale d’une ligne et d’un plan avec s = 12 et q = 1, 2, 3, 6.

Fig. 2 – Local 3D fractal dimension of a line and a plan with s = 12 and q = 1, 2, 3, 6.Fig. 2 – Dimension fractale locale d’une ligne et d’un plan avec s = 12 et q = 1, 2, 3, 6.

15The volume inside the cube depends on the coefficient h. The transformation into voxels can smooth the surface or, on the contrary, can enhance the surface roughness. When h = 1, the volume depends directly on the resolution of the DEM. Varying the coefficient h in Eq. (4) implies a change in the number of voxels vs, and thus the volume occupied inside the cube (fig. 1). High h values accentuate the surface roughness because the weight of the first term in this formula is more important than the second one (s/2). On the contrary, low values produce smooth surfaces because the second term of the formula is predominant. It is then possible to adapt the h value according to the nature of the DEM surface and the level of information expected.

16Determining the most appropriate range of grid sizes is a common problem for fractal dimension estimation (Foroutan-pour et al., 1999). The size of the frame used for the computation can be an even or an odd number; the range of grid size can be a power of two (Biswas et al., 1998) or a succession of integers. In our treatments, based on various tests applied to training zones, the size q of the sliding cube is chosen to be 12 or 24. Using these sizes, the slope P=ln(q)/ln(Ns) is more regular when exact dividers are employed. For a flat surface and a volume, the slopes P obtained are -2 and -3 respectively. Even if these proportions are respected and do not produce a different result after normalisation, odd sliding windows provide a slight deviation of the slope P for theoretical forms. The exact dividers, which define the range of grid sizes q, never correspond to half of the window side and then only a portion of the slope is taken into account in calculating the fractal dimension. When applying a 12 x 12 x 12 cubic pattern, five exact dividers (1, 2, 3, 4, and 6) can be found, excluding any border effect. In the second case (24 x 24 x 24), one can find seven exact dividers (1, 2, 3, 4, 6, 8, and 12). Treatments that employ these cubic sizes show that using the dividers 1, 2, 3, and 6 in the first case and the dividers 1, 2, 3, 6, and 12 in the second produces no bias. As the second term in Eq. 4 depends on the window size s, the use of the smaller window displays small irregularities; the larger window exhibits the general features of the DEM studied.

17The fractal dimension is calculated locally for each surface pixel in a sliding window centred on this pixel. As discussed above, the quality of results is greater with an even size than with an odd one. However, in the first case, the concept of ‘centre’ must be defined because the studied point cannot be located exactly in the centre but in one of the four points surrounding this centre. Using one of these points leads to the deviation of the result that favours one orientation.

18These remarks imply that the procedure must either calculate the local fractal dimension on each of these four pixels with the risk of blurring the response when calculating the mean, or define the value of the centre by taking the average value of the four pixels located at the centre of the window. Ic in Eq. (4) is calculated according to the second alternative in order to avoid blurring and to minimize the computation cost.

Test areas and results

19In order to test its ability to extract morphologic features according to different lithological and tectonic conditions, the local fractal dimension was applied in two regions. The first example, located in the region of Vittel (NE France), characterizes a sedimentary area affected by faults. The second example, located in eastern Anatolia (Turkey), is a volcanic region. In order to illustrate the performance of the procedure, the treatments concerning the Vittel region show mainly the variation of coefficient h and range of grid sizes effect, whereas the variation of the window size has been employed in the second example.

The Vittel region

Fig. 3 – Vittel region.
Fig. 3 – Région de Vittel.

Fig. 3 – Vittel region. Fig. 3 – Région de Vittel.

A: shaded relief map; B: Vittel geological sketch map: Symbols; B = Buntsandstein; M = Muschelkalk; K = Keuper; R = Rhetian; C: local fractal with s = 12, h = 100; D: local fractal with s = 12 and h = 2.
A: MNT ombré; B: carte géologique de Vittel; symboles: B = Buntsandstein; M = Muschelkalk; K = Keuper; R= Rhétien; C: dimension fractale locale avec s = 12 et h = 100; D: dimension fractale locale avec s = 12 et  h = 2.

20The Vittel area is located in northeastern France, 105 km west of the Rhine graben, between 48° and 49° N and 5° and 6° E. In the Vosges massif, outcrops of homoclinal Triassic and Liassic units overlie Variscan metamorphic and intrusive rocks. The east-west-trending Vittel fault is a major composite fault zone of complex geometry, resulting from the reactivation of a Palaeozoic discontinuity. In the region located between this fault and the N45-trending Esley fault, different groups of faults favour the water circulation that supplies the commercial Vittel spring water (Sykioti, 1994).

21The DEM of this area is an IGN (Institut Géographique National) product with a horizontal resolution of 50 m and a conical projection (fig. 3A). The corresponding geological map (Fig. 3B) has been geometrically corrected according to this projection in order to be overlaid on the DEM and compared with the local fractal results. The zone size is about 13.75 x 23.55 km (275 lines x 471 columns). Several treatments were applied to the DEM using a cubic size of 12 pixels and modifying only the value of the coefficient h. According to Eq. (4), decreasing coefficients contribute to smooth the roughness because a voxel inside the cube corresponds to a bigger hypsometric interval (ps/h). Then, the procedure detects only high altitude variations. On the contrary, when coefficient h increases, the surface roughness is accentuated and therefore small altitude variations are detected as well as large altitude ones.

22Thus, the position and extent of the Vittel fault (fig. 3C) are precisely highlighted and its prolongation westwards remains clear. When using high h coefficients (i.e., 100), every elevation difference is taken into account. Then, a thalweg feature in the testing cube corresponds to a filled cube with a gully whatever the altitude variation. The fractal dimension obtained is close to 3. On the other hand, a crest line is represented by a “wall” as formerly described (see Fig. 2) and its fractal dimension is close to 2 (fig. 3C). These results are obtained using a window size s = 12 and a range of grid sizes varying between 1 to s. These ensure the continuity of the features by producing a standardization of the estimated fractal dimension and increase the detection of crest and thalweg features.

23In contrast, the major units encountered in the studied zone were underlined using a low h coefficient but greater than one, because the region is relatively flat. With a coefficient equal to 2 (fig. 3D), the main structural and geological units are detected. A high value is obtained when, at the observation scale, the studied zone inside the moving window is relatively flat. The escarpment that corresponds to the Rhetian border is strong enough to induce in the moving window high roughness values codified with low grey tone values. Thus, the limit of the Rhetian appears clearly and corresponds exactly to the limits reported in the geological map. Concerning the Muschelkalk (limestones about 50 m in thickness) and the Bundsandstein (Vosgian sandstones), the incision depth of the drainage network is closely related to the respective nature of these two formations. For this reason, these units are clearly individualized and the contact corresponding to the N45-trending Esley fault is obvious, as well as the stratigraphical contact located at the base of Muschelkalk escarpment. The results are improved by superposing the geological map and the image resulting from the local fractal treatments. The frequency of the local fractal dimension depends on the sliding window size (fig. 4). Using a great size, the histogram is unimodal and becomes multi-modal with decreasing sizes, allowing us to obtain various classes corresponding to each mode.

Fig. 4 – Distribution of fractal dimension according to sliding window size.
Fig. 4 – Distribution de la dimension fractale en fonction de la taille de la fenêtre mobile.

Fig. 4 – Distribution of fractal dimension according to sliding window size.Fig. 4 – Distribution de la dimension fractale en fonction de la taille de la fenêtre mobile.

The Ararat volcano

Fig. 5 – Case study of Mount Ararat.
Fig. 5 – Exemple du Mont Ararat.

Fig. 5 – Case study of Mount Ararat. Fig. 5 – Exemple du Mont Ararat.

A) shaded relief map; B) geological map (after Yilmaz et al. 1998). 1: basalt, 2: hypersthene basalt 3: hyalobasalt, 4: andesite and associated pyroclastic rocks, 5: hypersthene andesite, 6: hyaloandesite, 7: moraines, 8: alluvial and glacial fans, 9: alluvium, 10: basement. I: Permanent ice cap; Q: Quaternary; C) local fractal with s = 12; D) local fractal dimension with s = 24.
A) MNT ombré; D) carte géologique (Yilmaz et al. 1998). 1: basalte; 2: basalte à hypersthène; 3: hyalobasalte; 4: andésite et roches pyroclastiques associées; 5: andésite à hypersthène; 6: hyaloandésite, 7: moraines; 8: cône de déjection; 9: alluvion; 10: complexe de base; I. neiges permanentes, Q:Quaternaire; C: dimension fractale locale avec s = 12; D: dimension fractale locale avec s = 24.

24With an elevation of 5123 m a.s.l., Mount Ararat is the largest volcanic centre and the highest point of eastern Anatolia. It corresponds to a compound stratovolcano formed by the ‘Greater’ and the ‘Lesser’ Ararat. This N150-trending, elongated massif lies inside a pull-apart basin with a similar trend (Pearce et al., 1990). This fault system is a horsetail splay structure and in addition to the right-lateral slip component, these faults have a normal component. The horsetail splay fault system cuts through the Greater and Lesser Ararat volcanoes and controls the position of the main eruptive centres of Ararat, as well as the alignment of parasitic cones (Karakharian et al., 2002). The volcanic activity of this region continued without interruption until historical times, possibly reaching its climax during the late Miocene and Pliocene (6 to 3 Ma). During Quaternary times, the volcanism appears to have been restricted to a few localised centres (Yilmaz et al., 1998). The youngest eruptions from the Ararat volcano are older than 10 000 years, and isotopic dating of the youngest lavas yields an age of 20 000 years (Yilmaz et al., 1998). Recently, geomorphic criteria applied to Mount Ararat have been used in order to study the tectonics of eastern Anatolia (Adiyaman et al., 2003). These criteria refer to the morphometric method developed by F. Garcia-Zuniga and J.-F. Parrot (1998) in order to define the base line and measure its elongation, and to the “Voxel wall” procedure (Baudemont and Parrot, 2000). The latter consists in parameterizing the DEM surface in a real 3D space using a vertical wall of voxels.

25The DEM of the Ararat zone (fig. 5A) is generated using the method developed by H. Taud et al. (1999). The pixel size equals 30 m. The size of this area is about 9.57 x 9.42 km (319 lines x 314 columns). The surface roughness of the different lava flows has been studied by means of different window sizes. A small window size (s = 12 with h = 4) reveals that the roughest surfaces mainly correspond to the more recent features, i.e., to the extent of the lava flows erupted from the recent parasitic cones located on the western flank. The extent of the eastern gully coming from the Greater Ararat is clearly observed, as well as the total base line of the whole edifice recently described by O. Adiyaman et al. (2003). A greater window size (s = 24 with h = 4) takes into account the regional structural features. As the local fractal values decrease in relation to brutal hypsometric variation as well as the character of the studied shape, a great window size reveals for instance the presence of isolated volcanic cones (low values), the presence of domes (high values surrounded by a crown of low values). The local fractal results depend on the size of the moving window as the frequency of the 2D signal is related to the surface roughness characterizing the type of studied lava flow. On the other hand, the value of the general roughness wavelength is an indicator of the regional tectonics. As illustrated by figures 5C and 5D, the SW quarter of the stratovolcano presents a globally lower roughness signature, the eastern limit of which emphasizes the presence of a large fault zone passing through the Greater Ararat crater. This fault corresponds to a branch of the horsetail splay structure described by A. Karakharian et al. (2002). In addition to the right-lateral slip, these authors assume that these faults have a normal component. The SW zone that has a weaker roughness response could indicate a huge collapse structure. In fact, the volcanic material (hypersthene andesite) observed in this zone (fig. 5B) is globally different from the material outcropping in the eastern part of the stratovolcano (mainly basalts and hypersthene basalts).

26Moreover, in order to study the effect induced by different window sizes and different values of h, the units of the geological sketch map drawn by Y. Yilmaz et al. (1998) have been reordered according to their petrographic nature. Six volcanic items (types of lava) are present: andesite, hyaloandesite, hypersthene andesite, basalt, hyalobasalts and hypersthene basalt. On the other hand, alluvium, fluvial and glacial fans, moraines and basement have been also taken into account in the statistical study. This geological map has been geometrically corrected in order to be overlain on the DEM and compared with the local fractal results. The mean of the local fractal dimension measured on the surface of each item is reported in the figure 6. The value of the fractal dimension depends on these two coefficients. Some items present specific values that allow us to characterize or to group them in different classes. When the window size is equal to 6 (fig. 6A), all the items except the alluvial formations are comprised in an unique class while the h value is equal to 1. For h values greater than 1, three classes can be observed: low fractal values are related to alluvium, basement and fans; median values to basalt, hyalobasalt and hypersthene andesite (the two formers are strongly connected); and high values to the remaining items. The same distribution is observed using window sizes of 24 where different values of h do not generate any variation (fig. 6C). On the contrary, with a 12 window size (fig. 6B) the items are classified differently. For instance, andesites are characterized by low values when s = 12 and high values when s equals 6 or 24. Therefore, using the mean value, statistical classifications can be obtained by taking into account the results produced by different values of coefficients s and h.

Fig. 6 – Evolution of mean fractal value of the studied items as a function of window size and coefficient h.
Fig. 6 – Evolutions de la moyenne de la dimension fractale des thèmes étudiés en fonction de la taille de la fenêtre et du coefficient h.

Fig. 6 – Evolution of mean fractal value of the studied items as a function of window size and coefficient h.Fig. 6 – Evolutions de la moyenne de la dimension fractale des thèmes étudiés en fonction de la taille de la fenêtre et du coefficient h.

27The items corresponding to large zones seem to be formed by a heterogeneous ensemble corresponding to different erupted materials and the presence of different volcanic flows. High values are registered at the bulged top of the lava flows because this top appears as a flat surface at the observed scale, while the lateral flanks present lower values. This is the case for the isolated flows such as the hypersthene basalt from the Lower Ararat (fig. 5C). Hence, the high roughness alignments detected by the treatment, which are not present in the different geological units mapped by Y. Yilmaz et al. (1998), could be related to different volcanic flows running from the eruptive centres.

Conclusion

28This paper aims at examining the relation between the surface roughness and the geological and geomorphologic features. It is well known that the fractal dimension can be used as a parameter characterizing the surface roughness and the landscape shape. The treatment proposed here measures locally this dimension by an adaptive box counting 3-dimensional method. The algorithm is based on the use of two parameters: the size s of the moving window and a scaling factor h defining the hypsometric interval taken into account. This treatment has been applied to a DEM surface of two regions characterized by a difference concerning their lithological and structural conditions. The first example is related to the study of a sedimentary and faulted region in the Vittel area. In this case, two main results have been obtained according to the different configurations of the parameters involved in the calculation. As a first result, the limits of the geological units detected by the local fractal dimension correspond to these limits as they have been mapped on the field. The position of these limits match either to the border of a homogeneous morphologic feature or to the transition zone between two different incision depths of the drainage network. On the other hand, by using low hypsometric intervals by means of the coefficient h, the procedure can detect fine features. We can therefore underline the Vittel fault as well as its western prolongation. The second example focuses on the volcanic region of Mount Ararat. In this case, the relation between the local fractal measurements and the different volcanic formations which have built up this stratovolcano has been studied by using a statistical approach. Different volcanic classes can be distinguished by using their mean fractal dimension value. On the other hand, the results provided directly by the local fractal dimension show the extent of the volcanic flows according to their surface roughness. Furthermore, a collapse structure related to the regional strike slip faulting has been detected. The results obtained with the two training DEMs corresponding to lithological conditions that are completely different, corroborate the efficiency of the procedure. The latter allows the retrieval of information about small structural features in flat zones as well as general morphologic characteristics of rugged areas.

29Based on the first results, one can conclude that the roughness of a surface is strongly correlated with the nature of the material that forms the geology of a studied region. In the case of the Vittel region mainly formed by a sedimentary sequence and studied as a training set, the application of the treatment by means of a descriptive analysis shows that this technique is efficient to extract and precise the meaning of the different morphometric features. In the case of more complex region such as the Mount Ararat that presents unknown structures, the extraction and the meaning of the morphometric features resulting from the treatment, needs to understand the nature of the filtering result; with such a goal, it is necessary to analyse precisely the results provided by using different values of the coefficients s and h, in order to realize an objective synthesis.

30Actually, the local fractal dimension provides useful information about geological and geomorphologic features. The proposed method detects different types of structures according to the observation scale and provides useful information concerning the geological and geomorphologic mapping.

Top of page

Bibliography

Adiyaman O., Parrot J.-F., Chorowicz J., Baudemont F., Köse O. (2003) – Geomorphic criteria for volcanoes from numerical analysis of DEMs. Application to the tectonics of Eastern Anatolia. Zeitschrift für Geomorphologie, 47, 235-250.

Baudemont F., Parrot J.-F. (2000) – Structural analysis of DEMs by intersection of surface’s normals in a 3D accumulator space. IEEE Transactions on Geoscience and Remote Sensing, 38, 1191-1198.

Biswas M. K., Ghose T., Guha S., Biswas P. K. (1998) – Fractal dimension estimation for texture images: A parallel approach. Pattern Recognition Letters, 19, 309-313.

Chase C. G. (1992) – Fluvial landsculpting and the fractal dimension of topography. Geomorphology, 5, 39-57.

Chaudhuri B. B., Sarkar N. (1995) – Texture segmentation using fractal dimension. IEEE Transactions on Pattern Analysis and Machine Intelligence, 17, 72-77.

Chen C. C., Daponte J. S., Fox M. D. (1989) – Fractal feature analysis and classification in medical imaging. IEEE Transactions on Medical Imaging, 8, 133-142.

Cheng Y. C., Lee P. J., Lee T. Y. (1999) – Self-similarity dimensions of the Taiwan Island landscape. Computers and Geosciences, 25, 1043-1050.

Clarke K. C. (1986) – Computation of the fractal dimension of topographic surfaces using the triangular prism surface area method. Computers & Geosciences, 12, 713-722.

Datcu M., Luca D., Seidel K. (1996) – Wavelet-Based Digital Elevation Model Analysis. 16th EARSeL (European Association of Remote Sensing Laboratories Symposium), Rotterdam, Brookfield, 283-290.

Dellepiane S., Giusto D. D., Serpico S. B., Vernazza G. (1991) – SAR image recognition by integration of intensity and textural information. International Journal of Remote Sensing, 12, 1915-1932.

Falconer K. (1990) – Fractal Geometry Mathematical Foundations and Applications. J. Wiley, Chichester, 288 p.

Foroutan-pour K., Dutilleul P., Smith D. L. (1999) – Advances in the implementation of the box-counting method of fractal dimension estimation. Applied Mathematics and Computation, 105, 195-210.

Garcia-Zuñiga F., Parrot J.-F. (1998) – Analyse tomomorphométrique d’un édifice volcanique récent: Misti (Pérou). Comptes Rendus de l’Académie des Sciences, 327, 457-462.

Huang J., Turcotte D. L. (1989) – Fractal mapping of digitized images: Application to the topography of Arizona and comparisons with synthetic images. Journal of Geophysical Research, 94, 7491-7495.

Karakhanian A., Djrbashian R., Trifonov V., Philip H., Arakelian S., Avagian A. (2002) – Holocene-historical volcanism and active faults as natural risk factors for Armenia and adjacent countries. Journal of Volcanology and Geothermal Research, 113, 319-344.

Klinkenberg B., Goodchild M. F. (1992) – The fractal proprieties of topography: A comparison of methods. Earth Surface Processes and Landforms, 17, 217-234.

Mandelbrot B. (1982) – The fractal geometry of nature. Freeman, San Francisco, 460 p.

Mark D. M., Aronson P. B. (1984) – Scale dependent fractal dimensions of topographic surfaces: An empirical investigation, with applications in geomorphology and computer mapping. Mathematical Geology, 16, 671-683.

Pearce J. A., Bender J. F., De Long S. E., Kidd W. S. F., Low P. J., Güner Y., Saraglu F., Yimaz Y., Moorbath S., Mitchell J. G. (1990) – Genesis of collision volcanism in Eastern Anatolia, Turkey. Special Publication, Journal of Volcanology and Geothermal Research, 44, 189-229.

Peitgen H. O., Jurgens H., Saupe D. (1992) – Chaos and Fractal: New Frontiers of Science. Springer, New York, 984 p.

Pentland A. (1984) – Fractal-based description of natural scenes. IEEE Transactions on Pattern Analysis and Machine Intelligence, 6, 661-674.

Polidori L., Chorowicz J., Guillande R. (1991) – Description of terrain as a fractal surface, and application to Digital Elevation Model quality assessment. Photogrammetric Engineering & Remote Sensing, 57, 1329-1332.

Russ J. C. (1999) – The image Processing Handbook. 3rd edition, CRC Press, Boca Raton, 771 p.

Sarkar N., Chaudhuri B. B. (1994) – An efficient differential box-counting approach to compute fractal dimension of image. IEEE Transactions on Systems Man and Cybernetics, 24, 115-120.

Sykioti O. (1994) – Méthodologie et imagerie numérique multisource de la surface topographique. Application à différents contextes hydrogéologiques: Vittel (Lorraine) et Verneuil-sur-Avre (Perche). Thèse de l’université Pierre et Marie Curie (Paris 6), 108 p.

Taud H., Parrot J.-F., Alvarez R. (1999) – DEM generation by contour line dilation. Computers & Geosciences, 25, 775-783.

Toennies K. D., Schnabel J. A. (1994) – Edge detection using the local fractal dimension. Proceedings IEEE Seventh Symposium, Computer-Based Medical Systems, 34-39.

Wilson J. P., Gallant J. C. (2000) – Terrain analysis. Principles and Applications. John Willey & Sons, New York, 479 p.

Top of page

Annex

Version abrégée

La rugosité ou la texture des Modèles Numériques de Terrain (MNT) est susceptible de fournir des informations relatives à la géologie régionale. En effet, les MNT étant une représentation de la surface, différents attributs peuvent les décrire. Entre autres, la dimension fractale permet de caractériser la texture dans la mesure où la topographie terrestre est sensée présenter un comportement fractal indépendamment de l’échelle d’observation. Cet article concerne l’étude de la rugosité de surface à l’aide de la mesure de la dimension fractale locale dans un espace tridimensionnel, en vue de mettre en évidence ou d’accentuer divers traits géomorphologiques.

La géométrie fractale est une description mathématique des formes naturelles. Un objet fractal est trop complexe pour être décrit dans un espace cartésien. Seule la dimension fractale est à même de mesurer un objet complexe. Différentes méthodes ont été proposées pour mesurer cette dimension. L’une d’entre elles est largement utilisée : il s’agit du « comptage de boîtes » pouvant être appliqué à tous type de formes, fractales ou non.

Dans le cas présent, nous proposons de mesurer localement au sein d’un cube la dimension fractale en se fondant sur cette méthode. L’avantage de ce traitement réside dans le fait que les « voxels » décrivant le volume pris en compte pour faire ce calcul dépendent directement de l’altitude des pixels décrivant la surface du MNT. Le « voxel » est lui-même un cube dont la base est un pixel et la hauteur une tranche d’altitude correspondant à la dimension du coté du pixel. Une altitude donnée est elle ainsi représentée par un empilement de voxels ou cubes élémentaires. En fait la variation du coefficient h décrit plus loin permet de modifier la hauteur de cette tranche d’altitude. La procédure est la suivante : au sein d’un cube de taille s x s x s centré sur un pixel décrivant la surface du MNT, le volume correspondant à la section du MNT prise en compte est un ensemble de voxels dont le nombre est compris entre 0 et s. Le nombre de voxels présents dans le cube est égal à :

vs peut se calculer facilement de la manière suivante en se plaçant dans l’espace bidimensionnel :

et  

I correspond à l’image traitée, Ic est la valeur du pixel central, ps la taille du pixel, h un coefficient définissant la résolution verticale et s la taille du cube. Ce cube est divisé en boîtes dont la taille varie de 1 à s/2. Chaque boîte est considérée comme remplie quand elle contient au moins un voxel. En d’autres termes, la valeur maximale de vs(i,j)/q est calculée en appliquant l’équation suivante :

La dimension fractale correspond à l’inverse de la pente P=ln(q)/ln(Ns), q étant la taille de la boîte et Ns le nombre total de boîtes remplies.

Les traitements effectués à titre d’exemple concernent deux milieux choisis pour leurs différences tectoniques et lithologiques et s’appliquent à des zones test antérieurement étudiées, en vue de valider la méthode. Ils illustrent quelques-unes des possibilités qu’offre cette approche au plan morphologique et structural. Le premier, situé dans la région de Vittel (France), est une zone sédimentaire faillée. Le second, situé en Anatolie orientale (Turquie), concerne l’étude d’un strato-volcan complexe, le Mont Ararat. Dans les deux cas, des études comparative et statistique prenant en compte les données existantes et les résultats issus du filtrage provenant de la mesure de la dimension fractale locale, ont été réalisées en vue d’évaluer les résultats que fournit cette approche.

Dans le premier cas, la relation existant entre la géomorphologie et les résultats obtenus à l’aide du filtre décrit antérieurement est particulièrement nette. Il convient de noter que l’on est ici en présence d’une série monoclinale, la surface de chacune des couches géologiques répondant en fonction de ses caractéristiques propres. Les différentes unités stratigraphiques sont mises en évidence et cernées, soit par le biais des épaulements qui les limitent comme cela peut s’observer dans le cas du Rhétien, soit par l’importance et la profondeur des incisions qu’engendre le réseau hydrographique dans ces formations. Dans ce dernier cas, les formations du Muschelkalk et du Buntsandstein sont clairement définies en utilisant ce critère, ainsi que l’accident N-S correspondant à la faille d’Esley qui les met en contact. Par ailleurs, en jouant sur les valeurs du coefficient h, il est possible de détecter les lignes de crête et les thalwegs et ainsi de souligner la faille de Vittel et son prolongement occidental.

L’utilisation de la dimension fractale locale dans le cas du Mont Ararat est d’un usage plus délicat dans la mesure où un strato-volcan est une structure complexe, tant au plan de la nature du matériel volcanique que de la tectonique. En fait, la dimension fractale locale calculée en utilisant des fenêtres de grande taille se révèle utile pour mettre en évidence les grandes unités structurales qui caractérisent ce massif, ainsi que l’extension des différentes coulées volcaniques qui le constituent. Situé dans un réseau de failles de direction N-S dont l’une passe par le Grand Ararat, l’ensemble présente en premier lieu une zone d’effondrement SW induisant une réponse de la rugosité de surface différente de celle qui domine sur le reste du massif. Cette différence se retrouve au niveau de la lithologie des émissions effusives. Un filtrage de grande taille révèle et précise la nature des édifices (localisation de la ligne de base, position des cônes de cendres et/ou des dômes volcaniques). De plus, les coulées répondent comme des linéaments caractérisés par une haute valeur au toit de la formation (absence de rugosité au sein de la fenêtre d’observation), cernés par des valeurs plus faibles correspondant aux flancs latéraux, en raison de l’écart hypsométrique enregistré qui abaisse la valeur de la dimension fractale. Il est ainsi possible de définir au sein des grands ensembles pétrographiques une partie des éléments qui les composent.

Le résultat de l’application de ce traitement sur ces deux zones différentes nous conduit aux conclusions suivantes. Lorsque la zone d’étude correspond à un ensemble tectoniquement homogène, la compréhension de la réponse passe par une analyse descriptive relativement facile à réaliser et le traitement peut donc représenter dans ce cas un outil efficace permettant de préciser les traits géomorphologiques. En revanche, dans le cas d’ensembles complexes, la signification des traits morphologiques mis en évidence par cette approche méthodologique nécessite une étude de détail de la réponse de tous les objets mis en jeu, dans le but de réaliser une synthèse objective.

Top of page

List of illustrations

URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-1.png
File image/png, 4.7k
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-2.png
File image/png, 5.3k
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-3.png
File image/png, 4.2k
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-4.png
File image/png, 8.4k
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-5.png
File image/png, 5.5k
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-6.png
File image/png, 9.1k
Title Fig. 1 – Local fractal calculation by using a 3-dimensional box counting. Fig. 1 – Calcul de la dimension fractale locale par la technique tridimensionnelle du « comptage de boîtes ».
Caption Example with s= 12 and different grid size q. A: initial topography equivalent to q= 1; B: q= 2; C: q= 3; D: q= 6.Exemple avec s = 12 et différentes tailles de maille q. A: topographie initiale équivalente à q = 1; B: q = 2; C: q = 3; D: q = 6.
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-7.png
File image/png, 50k
Title Fig. 2 – Local 3D fractal dimension of a line and a plan with s = 12 and q = 1, 2, 3, 6.Fig. 2 – Dimension fractale locale d’une ligne et d’un plan avec s = 12 et q = 1, 2, 3, 6.
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-8.png
File image/png, 44k
Title Fig. 3 – Vittel region. Fig. 3 – Région de Vittel.
Caption A: shaded relief map; B: Vittel geological sketch map: Symbols; B = Buntsandstein; M = Muschelkalk; K = Keuper; R = Rhetian; C: local fractal with s = 12, h = 100; D: local fractal with s = 12 and h = 2.A: MNT ombré; B: carte géologique de Vittel; symboles: B = Buntsandstein; M = Muschelkalk; K = Keuper; R= Rhétien; C: dimension fractale locale avec s = 12 et h = 100; D: dimension fractale locale avec s = 12 et  h = 2.
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-9.png
File image/png, 400k
Title Fig. 4 – Distribution of fractal dimension according to sliding window size.Fig. 4 – Distribution de la dimension fractale en fonction de la taille de la fenêtre mobile.
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-10.png
File image/png, 4.9k
Title Fig. 5 – Case study of Mount Ararat. Fig. 5 – Exemple du Mont Ararat.
Caption A) shaded relief map; B) geological map (after Yilmaz et al. 1998). 1: basalt, 2: hypersthene basalt 3: hyalobasalt, 4: andesite and associated pyroclastic rocks, 5: hypersthene andesite, 6: hyaloandesite, 7: moraines, 8: alluvial and glacial fans, 9: alluvium, 10: basement. I: Permanent ice cap; Q: Quaternary; C) local fractal with s = 12; D) local fractal dimension with s = 24.A) MNT ombré; D) carte géologique (Yilmaz et al. 1998). 1: basalte; 2: basalte à hypersthène; 3: hyalobasalte; 4: andésite et roches pyroclastiques associées; 5: andésite à hypersthène; 6: hyaloandésite, 7: moraines; 8: cône de déjection; 9: alluvion; 10: complexe de base; I. neiges permanentes, Q:Quaternaire; C: dimension fractale locale avec s = 12; D: dimension fractale locale avec s = 24.
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-11.png
File image/png, 184k
Title Fig. 6 – Evolution of mean fractal value of the studied items as a function of window size and coefficient h.Fig. 6 – Evolutions de la moyenne de la dimension fractale des thèmes étudiés en fonction de la taille de la fenêtre et du coefficient h.
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-12.png
File image/png, 11k
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-13.png
File image/png, 4.0k
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-14.png
File image/png, 6.7k
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-15.png
File image/png, 4.0k
URL http://journals.openedition.org/geomorphologie/docannexe/image/622/img-16.png
File image/png, 6.9k
Top of page

References

Bibliographical reference

Hind Taud and Jean-François Parrot, “Measurement of DEM roughness using the local fractal dimension”Géomorphologie : relief, processus, environnement, vol. 11 - n° 4 | 2005, 327-338.

Electronic reference

Hind Taud and Jean-François Parrot, “Measurement of DEM roughness using the local fractal dimension”Géomorphologie : relief, processus, environnement [Online], vol. 11 - n° 4 | 2005, Online since 01 January 2008, connection on 29 March 2024. URL: http://journals.openedition.org/geomorphologie/622; DOI: https://doi.org/10.4000/geomorphologie.622

Top of page

About the authors

Hind Taud

Jean-François Parrot

By this author

Top of page

Copyright

The text and other elements (illustrations, imported files) are “All rights reserved”, unless otherwise stated.

Top of page
Search OpenEdition Search

You will be redirected to OpenEdition Search