1Volcanic cinder cones formed by cinder and pyroclastic debris are generally considered as truncated cones with a crater located in the summit (Macdonald, 1972). The first geomorphological studies have shown that morphologic changes occur with time and are able to provide information about the age of the edifice (Colton, 1967; Scott and Trask, 1971). Porter (1972) was the first to define quantitative ratios between different parameters in order to characterize the volcanic shape: i.e., the ratio height of the cone versus the base diameter would be equal to 0.18 and the ratio between the crater diameter and the base diameter would remain at 0.40. Bloomfield (1975), using radiometric age determinations, observed that the first ratio decreases from 0.21 until 0.10 with time, meanwhile the second one increases from 0.40 to 0.83. On the other hand, according to Settle (1979), the shape characteristics of the volcanic cones are related to the nature of the material involved in the effusive process, and to the nature and duration of the erosion activity. Wood (1980a, 1980b) confirms and formalizes the morphometric parameters proposed by Porter (1972). Until now, numerous geomorphological studies concerning the geomorphic definition of the volcanic characteristics or the effect of the erosion processes, are based on these parameters and ratio (Dohrenwed et al., 1986; McFadden et al., 1986; Hasenaka, 1994; Noyola-Medrano et al., 1994; Hooper, 1995; Luhr et al., 1995; Hooper and Sheridan, 1998; Rech et al., 2001; Aranda-Gomez et al., 2003; Nemeth et al., 2005), even if their use remains problematic when the studied cone does not present a crater (Hasenaka and Carmichael, 1985a).
2The possibility to obtain more detailed geomorphic information has been explored by Garcia-Zuniga and Parrot (1998) who proposed using a Digital Elevation Model (DEM) and to define pattern recognition parameters applied to hypsometric slices describing the volcanic cone, from its base line to its summit. This approach described as tomomorphometric analysis registers the morphologic changes taking into account parameters such as the convexity index, direction of the principal axis, etc. This recent approach has been used to study the lithospheric motion of the Somalian and Arabian plates (Collet et al., 2000), the Anatolian volcanic massif (Ozlem et al., 2003) and the Chichinautzin volcanic cinder cone field, Mexico (Noyola and Parrot, 2005).
3Square-grid digital elevation models (DEMs) represent important and accurate tools to underline the different regional geomorphic features and to simulate various scenarios, as the possibility of storage and advances in computing technology increased strongly in recent years. The horizontal and vertical resolutions are sufficient to accurately calculate different parameters extracted from the DEM surface. The digital terrain analysis (Wilson and Gallant, 2000) allows defining primary and secondary attributes. The secondary attributes are devoted to estimate the role played by topography in the distribution of soil water or on the susceptibility of landscapes to erosion, for instance. Most of the primary attributes such as slope, aspect, plan and profile curvature are computed directly or by fitting an interpolation function to the DEM in order to calculate them (Moore et al., 1993b; Mitasova et al., 1996; Florinsky, 1998). These different attributes provide numerous geomorphological indicators (Tribe, 1991). They are used to describe the morphometry, catchment position, stream channels, etc., and to compute topographic attributes (Jenson and Domingue, 1988; Dikau, 1989; Moore et al., 1993a; Dymond et al., 1995; Giles, 1998; Borrough et al., 2000). An exhaustive survey of the publications concerning this topic can be found in the Pike’s report (2002). Various softwares such as TAPES-G (Moore, 1992) have been developed with such a goal.
4This paper aims at defining and computing parameters of volcanic cones in order to enable the measurement of evolution stages of these cones due to the erosion, catastrophic events or human activity. The first approach consisted of defining some characteristic features that can be useful for this purpose. The algorithm presented in this paper has been developed in C++. The automated parameterisation of volcanic cones is based on seven parameters: the volume of the volcanic cone, the volcanic base line radius and eventually its elongation, the total height of the cone from the base line until its summit, the crater radius when existing, the depth of this crater, the mean dipping angle outside the crater and the mean dipping angle inside the crater.
5All these parameters and their relationships are not only able to characterize a volcanic cone, but can be also used to reconstitute the original volcanic landform. Moreover, the comparison between the reconstituted cone and the presently observed shape allows us to assess the erosion and evolution stages. The developed algorithm is based on a calculation that only requires the altitude of the volcanic base line and the coordinates of a point considered as the center of the volcanic cone when it is present or to the summit if this feature is lacking. The algorithm is described in the following section. In the second section, the procedure that allows us to reconstitute the cone taking into account the formerly obtained parameters, is presented. Finally, the results are discussed in the third section.
6The main lines of the computation consists in defining and extracting the parameters that characterize the studied volcanic cone and, as a second step, in reconstructing the original landform of the studied volcanic cone. The difference observed between these two cones also corresponds to a parametrization of the erosion processes.
7The first calculated parameters are the volcanic height H of the cone, the base line radius BR which present two components: the minimum base line radius BRmin and the maximum base line radius BRmax allowing calculation of the elongation ratio ER of the volcanic base line, the depth of the crater HC as well as the crater radius CR and its two components: the minimum crater radius CRmin and the maximum crater radius CRmax.
8The second group of parameters concerns the mean dipping angle ∝ outside the crater and the mean dipping angle inside the crater. One can notice that contrary to the procedure proposed by Wood (op. cit.), these values result from the slope calculation in each DEM’s point. Complementary information is provided: the percentage of downward angle located on the cone flank, as well as the percentage of the upward angle inside the crater, in order to measure irregularities that can emphasize for instance the presence of eroded zones (natural or artificial) or to control the degree of the smoothness of the cone. The third generation of parameters concerns the calculation of the bi and tri-dimensional surface, as well as the volume of the edifice.
9At the beginning of the procedure, the user has to define the altitude of the base line and the coordinates of the volcanic center CE. As the first step consists of researching on the DEM all the pixels whose altitude is equal or superior to the altitude of this base line, a labellization is required in order to choose among the different pixel components, the surface corresponding to the first altitude slice of the studied edifice. This surface is the base on which the different parameters will be calculated.
10This measurement needs to research the center of mass CM of the base surface and the Principal Axis (PA) passing through this point. PA is calculated as follows:
tg(2∝) = 2μxy/(μyy - μxx) if (μyy - μxx) ≠ 0
11On the other hand, the center of mass CM (Xc, Yc) and the moments of second order xx , yy and xy are respectively equal to:
12Nbp is the number of pixels of the object and Xi Yi the coordinates of the pixel i.
The greater distance between CM and the crossing point between the perimeter of the surface and PA corresponds to the value of the radius BRmax. The normal of PA passing through CM is calculated in order to define the value of BRmin. The elongation is then calculated.
13The calculation of the crater radius CR takes into account the estimated center of the edifice CE (coordinates defined by the user) and the highest altitude point. The normal of the line linking these two points allows researching the perpendiular radius, taking into account the highest point located on the both sides of this normal. The second value is either greater or lower than the first one. CRmax and CRmin depend on this comparison and then the elongation of the crater is computed. One can notice that a difference between the coordinates of CE and CM previously defined represents another morphometric parameter, i.e. the existing tectonic constraints during the volcano formation.
14The crater depth corresponds to the hypsometric difference between the altitude of the point CE and the highest altitude of the volcano. The height of the cone is the difference between this latter point and the altitude of the base line.
15The procedure used to calculate the volume consists of computing the volume of each altitude slice (in meters or decimeters according to the hypsometric resolution of the DEM). Each slice corresponds to the pixels whose altitude A’ is equal or greater than the former slice of altitude A that plays the role of the floor. The difference (A’-A) corresponds to the hypsometric resolution of the DEM. When the calculation is done on a slice, the roof of this slice plays the role of the floor for the following upper slice. We have to notice that the volumetric excess corresponding to the upper half part of a slice is actually balanced by the loss registered in the lower half part, blurring the staircase effect inherent in this procedure.
16Two types of surface can be defined, as well as the difference they present. The first calculated surface is the topographic surface directly calculated on the surface defined by the base line. In a first approximation (Pratt, 1978), this surface S2 is equal to:
17S2 = ΣPs + Σ (Pp / 2)
18where Ps and Pp are respectively the pixels describing the surface and the pixels belonging to the perimeter (fig. 1). This measure corresponds only to the surface occupied on the map, but it is also possible to measure the tri-dimensional surface S3.
Fig. 1 – Perimeter and surface pixels of an hypsometric slice.
Fig. 1 – Pixels du périmètre et de la surface dans une tranche altimétrique.
19In order to calculate the 3D surface, a new algorithm has been developed; each pixel is divided in eight rectangular triangles that converge in the pixel center. Figure 2 illustrates the computation process taking into account the altitude value of the studied pixel and the altitude values of the eight surrounding pixels. The altitude value of the pixel corners and the altitude of the center of each pixel side result from a linear interpolation between the altitude value of the pixel center and the altitude of the neighboring pixels.
Fig. 2 – Computation of the “3D surface” inside a pixel.
Fig. 2 – Calcul de la “surface tridimensionelle” à l’intérieur d’un pixel.
Fig. 3 – The 22 first templates, their corresponding bi-dimensional surface values inside the central pixel, and the corresponding triangular zones inside the central pixel that will be used in order to calculate the surface of an hypsometric slice.
Fig. 3 – Les 22 premiers patrons permettant de définir la portion d’un pixel du périmètre entrant dans le calcul de la surface d’une tranche d’altitude.
20The ratio S3D/S2 is another indication concerning the smoothing degrree of the studied cone. It becomes lower when the cinder cone is recent and does not present gullies, as well as when reconstituing the former shape (see table 1).
Table 1 – Parametric values of the Jocotitlán volcano (original form and reconstituted edifice).
Tableau 1 – Valeur des paramètres du volcan Jocotitlán (forme actuelle et édifice reconstitué).
* indicates that the involved algorithm encountered only one summit; for this reason, the two radius of the crater are similar.
* indique que la séquence algorithmique n’a rencontré qu’un seul sommet ; pour cette raison, les deux rayons du cratère ont la même valeur.
21The mean slope calculation is as follows. Following each line issued from CE (coordinates defined by the user) and linking a pixel belonging to the perimeter describing the base line, the algorithm compares the altitude of the successive points located on this line in order to know if the calculated angle corresponds to an ascending () or descending (∝) slope. A clockwise linking is done in order to scan the whole edifice and the local result is plotted if a previous one does not exist. On the other hand, the presence of local ascending slopes in the volcanic flank, as well as the presence of local descending slopes inside the crater are researched in order to control the regularity of the conic shape.
22In addition, the developed algorithm that can be provided by the author, reconstitutes the original volcanic cone using the following considerations. As formerly described a volcanic cinder cone corresponds to a truncated cone with a crater. The flank of the edifice is quite regular and smooth. Taking into account the value of the crater radius (CR) and the altitude of the highest point (HP), a circle is drawn whose altitude is equal to HP. If the coordinates of the center defined by the user are correct, this circle has to recover all the remnants of the crater wall.
23The calculated crater circle corresponds to the spatial reference on which the computation is based. Inside the crater, circular curve lines are drawn, the altitude of which are comprised between the altitude HP and the altitude value of the crater bottom; similarly on the flank of the volcano, circular curve lines are defined between the crater circle and the base line. The resulting DEM generation is obtained by using a curve dilation procedure (Taud et al., 1999).
24When the original altitude in some places is higher than the interpolated value, the original hypsometric value is preserved in order not to blur these features that can correspond to adventive effusive centers or local tectonic events. The generation of a reconstituted edifice allows the quantification of the differences registered by the volcanic cone due to erosion processes, collapses or human activity. The parameters, formerly described in the case of the original shape of the volcanic cone, are calculated in the same way in the case of the reconstituted cone. Their comparison is a key to underline the features resulting form erosion processes, collapses or human utilization.
25This paper was focused on two volcanic regions in Mexico. The first one comprises the Jocotitlan volcano, and the second is formed by two volcanoes from the Chichinautzin range.
26The Jocotitlán volcano (3950 meters) is a typical stratovolcano and dome complex, mostly composed by dacite lava flows. It is located in the central part of the Trans Mexican Volcanic Belt (TMVB), 60 km north of Mexico City. A huge collapse affected the NE sector of this edifice in pre-historic times. The associated avalanche covered an area of 80 km2 with a maximal runout distance of 12 km and an estimated volume of 2.8 km3 (Siebe et al., 1992). The presence of a major normal fault interesting the volcano suggests an extensional tectonic stress regime that could be responsible for triggering such a slope failure. The objective of this investigation was not directed towards studying the relationships between the volcano and the regional tectonic pattern, but to present and discuss the parameters calculation procedure (table 1), as well as the information obtained in terms of volume of the displaced material, and other measurements associated with the shape of the edifice. A 11-m resolution DEM produced taking into account a multidirectional interpolation (Parrot and Ochoa-Tejeda, 2004) method of curve lines was used as a base (fig. 4).
Fig. 4 – Contour lines of the Jocotitlán volcano used to produce the DEM.
Fig. 4 – Courbes de niveau du volcan Jocotitlán utilisées pour créer le MNT.
27By applying the proposed calculation procedure, it was possible to precisely characterize the volcanic original shape. They are compatible with the former parameters proposed by Porter (1972) and Wood (1980a, b). For instance, the ratio Hco/Wco (Height versus Diameter) is equal to 0.17, and the main dipping angles (17.47 for the flank and 18.52 for the crater) are comprised within the range defined by these authors for such a type of volcano. The calculation of the ratio Wcr/Wco (Crater diameter versus Base line diameter) depends on the position of the base line. In this case, the lower value obtained for this ratio (0.21) suggests that the lower limit chosen for the volcanic complex is located below the original base line in order to calculate the volume mobilized by the collapse event. The relation 3D surface versus the global volume is an indicator of the regularity of the edifice. Considering the original form, this ratio is equal to 3.77, whereas the value decreased to 3.40 after the volcano reconstitution (fig. 5 and fig. 6).
Fig. 5 – Shadowed DEM of the Jocotitlán volcano. A: original shape; B: reconstituted volcano.
Fig. 5 – MNT avec estompage du volcan Jocotitlán. A : forme actuelle ; B : forme reconstituée.
Fig. 6 – 3D diagram of the Jocotitlán volcano. A: original shape; B: reconstituted volcano.
Fig. 6 – Bloc diagramme du volcan Jocotitlán. A : forme actuelle ; B : forme reconstituée.
28Moreover, such parameters used to describe the volcanic edifice can be also utilized to calculate the volume of the displaced material at the time of the collapse. The volume of the displaced material was calculated as 1.976 km3, that means a smaller volume than the assessment proposed by Siebe et al. (1992). The later took into account the volume of the scattered debris-avalanche deposits expressed by a hummocky topography (see fig. 5). In addition to the analysis of the changes of the original shape and the reconstituted volcanic edifice, this procedure also can be used to define different zones affected by erosion and further calculate the local dissection depth (fig. 7).
Fig. 7 – Dissection depth provided by the altimetric difference between the original form and the reconstituted edifice.
Fig. 7 – Profondeur de dissection calculée en fonction de la différence hypsométrique entre la forme actuelle et la forme reconstituée.
1: 1-47 m; 2: 48-94 m; 3: 95-141 m; 4: 142-188 m; 5: 189-235 m.
1 : 1-47 m ; 2 : 48-94 m ; 3 : 95-141 m ; 4 : 142-188 m ; 5 : 189-235 m.
29On the other hand, two volcanoes of the Chichinautzin range were also studied. They are located 40 km south of Mexico City. The first one, El Tezoyo, a volcanic dome, is mainly composed by scorias and it lacks of a crater. The second one, Volcan del Aire, is a cinder cone with a well preserved circular crater. Same calculation procedure was applied for the two volcanoes by using a 10-m resolution DEM (tab 2).
Table 2 – Parametric values of El Tezoyo and Volcan del Aire volcanoes (original form and reconstituted edifice).
Tableau 2 – Valeur des paramètres du volcan El Tezoyo et du volcan del Aire (forme actuelle et édifice reconstitué).
30In order to enhance visual effects, resulting treatments of both were included within the same image (fig. 8a and fig. 8b).
Fig. 8 – 3D diagram of El Tezoyo and Volcan del Aire volcanoes. A: original shape; B: reconstituted volcanoes.
Fig. 8 – Bloc diagramme des volcans El Tezoyo et Volcan del Aire. A : forme actuelle ; B : forme reconstituée.
31The analysis was performed on the upper part of the two volcanos, taking into account the lower closed contour line as the base line. The treatment was done in order to illustrate the peculiarity of features. For instance (table 1), the lower and greater radius of the crater remain the same. Consequently, such result means that the involved algorithm encountered only one summit, as it is the case for the Jocotitlan volcano. On the contrary, when the crater is almost complete (Volcan del Aire), these values are different (see tab 2).
32Moreover, the value of the flanks mean slope can be regarded as an indicator of the volcanic morphology (14° for the cinder cone, 10° for a dome structure and 17° for a stratovolcano complex).
33The procedure presented in this paper calculates directly from a DEM, different parameters in order to characterize a volcano. This fast computation procedure (1) can be applied on a volcanic cone independently of the presence or the absence of a crater, (2) it does not only concern regular cinder cones, and (3) it can be used to study a composite volcano. The first selected test zone, the Jocotitlán volcano, corresponds to a volcanic structure of about 1200 meters height covering a surface of 62.5 km2. An important collapse of the northern flank occurred in pre-historic times. The values obtained by applying such procedures are useful to make calculations of parameters such those proposed by Wood (1980). However, the method also can be succesfully applied to perform measurements and analysis such as 3D surface or edifice volume.
34The reconstitution of the former volcanic landforms offers the possibility to accurately compute the “erosion” and the volume of the displaced material. The whole measurement set, can be considered as a new tool, able to discriminate different volcanic cones located in a volcanic field, to quantify their relative ages and emphasize the degree of degradation due to erosive processes.