Volcaniclastic Deposit Analysis Using Grain Size, Digital and Stereological Techniques ()
1. Introduction
The volcanoclastic (or volcaniclastic) term was introduced [1] to group all clastic rocks, with a volcanic origin, into a single system; subdividing them into three types: 1) pyroclastic, 2) autoclastic, and 3) epiclastic [1] [2] . It is common to use the term “volcanoclastic deposits” for those formed by mainly sedimentary processes during a volcanic event, and “pyroclastic deposits” for those derived directly from volcanic eruptions; but in recent work, the term “volcanoclastic” is used as a general classification; even in new classification attempts dividing them as primary and secondary to distinguish between direct volcanic eruption products, or as reworked ones, respectively [3] .
Research on volcanoclastic deposits is relatively young. Since its inception, it has been examined following traditional sedimentological procedures, with the same problems as in the study of rocks and sedimentary deposits, where the grain shape and size, the geometry of the deposit and its internal fabrics are used as tools to determine the physical processes controlling their formation and deposition [4] . Although the determination of these stratigraphic relationships is generally the first step in almost all the volcanoclastic materials’ studies, petrological and petrographic principles, as well as geochemical techniques, are equally applied [5] .
The difference between a primary and a secondary volcanoclastic deposit is found within their facies [6] . A deposit showing evidence of a pyroclastic fall accumulation is considered as primary, and is subdivided into pyroclastic fall, pyroclastic density current, peperite (a sedimentary rock with fragments of younger igneous rocks formed when magma gets in contact with wet sediments [7] , and hyaloclastite (a volcanoclastic breccia composed of fragments and glass that formed during a submarine eruption [6] [8] . Meanwhile, secondary deposits indicate gravitational action or removal processes and successive re-sedimentation, and are subdivided into lahars and debris avalanches or, in case of reworking by surface processes, into epiclastic deposits. Although some proposals for more subdivisions [3] can be found in the literature, these tend to be specific in their use.
The term “lahar” is a geological term originated in Indonesia, used to refer to water-borne volcanic breccias (or a “mudflow”), without necessarily being of volcanic origin [9] . The term has become synonymous with volcanic debris flow [5] . Although it has been widely proposed to drop the term entirely due to its ambiguity by using “debris flow” instead [10] , however, the word “lahar” has been retained. This research follows previous work [8] , adopting it as a general term, and using “debris flow” and “hyper-concentrated flow” as flow types within a lahar.
A volcanic debris avalanche (VDA) is known as the process formed due to the disaggregation of a landslide, or landslide from a volcanic building [11] [12] . Both VDAs and lahars are considered as “volcanic mass-wasting flows”, or as movements of large amounts of detrital material down the side of a volcanic building due to gravity [13] .
Based on the Udden-Wenworth scale [14] [15] used to describe the size of sedimentary fragments, several scales have been developed for volcanoclastic deposits. Although they differ in their boundaries, most classifications separate clasts (or volcanoclasts) into three main types: 1) ash (<2 mm), 2) lapilli (2 - 64 mm), and 3) bombs or blocks (>64 mm), following the terms established [1] . It is common, especially in size distribution studies, to use Krumbein’s phi (ϕ) notation [16] for those limits. Qualifiers such as coarse, medium, and fine, to each clast size were added later [8] .
Nevertheless, there is no similar classification for secondary volcanoclastic deposits, using instead sedimentary terms to describe them, and naming as “volcanoclast” for any rock fragment greater than 2 mm in diameter, and matrix for the finest material. In turn, a distinction was made between an interclast matrix, supporting the volcanoclasts, and an intraclast matrix, which is part of the volcano clasts [17] .
Fabrics in volcaniclastic rocks and deposits are determined almost exclusively by their depositional mechanisms, and by measuring their surface distribution [5] and contrasting them with observations in nature and experiments [18] . Quantitative analyses of sediment size distributions are thus necessary, for a detailed description of the sediments [19] .
For unconsolidated sediments, different measurement methods are used depending on their maximum size. For coarse sediments (gravel and sand particle size), sieving is usually the most applied technique. However, in deposits where clasts exceed the block size, this technique becomes impractical, while other techniques such as clast counting and adhesive sampling, increase both the analysis’ time and cost.
During the late 1990’s computer-based image analysis techniques, were already emerging as a faster and cheaper alternative for several types of geological studies [20] . Since then, different optical techniques for granulometric analysis have been applied to rocks and sedimentary and/or volcanoclastic deposits, obtaining results just as accurate as traditional techniques [21] [22] [23] [24] . Its use also causes less alteration to the studied surfaces (common during sieving), as well as solving issues caused by large blocks being present, or in sampling areas that are difficult to get to.
Most optical methods are based on measuring the surface area of features in 2D slices, segmenting the image(s), and automatically producing a binary image, measuring each object within it [24] . The next step is to convert these 2D measurements to 3D data, by applying stereology [20] , which is “the study of three-dimensional (3D) structures, and particularly the measurement of key parameters describing such structures based on in two-dimensional (2D) images” [25] .
The present study conducts a sedimentological and grain-size analysis on volcanoclastic deposits, employing computerized granulometric analysis techniques based on orthogonal images to ascertain grain size distribution. The aim is to propose a genetic classification for these deposits.
2. Location Site
The study area consists of two sites located in the SW part of Chihuahua city, Mexico, nearby a well-known shopping mall development (Figure 1). Almost the entirety of Chihuahua city is located within the physiographic province of
![]()
Figure 1. Study area sites, SW part of Chihuahua city, Mexico, where “site A” (PA) and “site B” (PB) are located.
Basin and Range, specifically on the Chihuahua-Sacramento Valley, composed mostly of rhyolites, rhyolitic tuffs, and Tertiary carbonates [26] . The urban area is developed over Neogene and Quaternary sediments and sedimentary rocks.
The Mexican Geological Service reports in its H13-C66 (scale 1:50,000) geological chart volcanic and pyroclastic rocks of rhyolitic to intermediate composition in the study areas [27] .
3. Methodology
The study of volcanoclastic rocks and deposits is essentially a sedimentological study because their classification was built based on sedimentary material analyses [1] [2] . This was recognized and a methodology for the analysis of modern volcanoclastic deposits [4] was conceived. They listed a series of commonly measured properties, such as: thickness, maximum grain size, grain size distribution, proportions of the components, crystal or pumice content, density, and porosity. These parameters allow to distinguish between the eruptive units of a deposit and facilitate their correlation.
Although this work takes the quantitative analysis methodology approach, a different one had to be employed during field work at each outcrop, according to the availability and accessibility. A crucial aspect was the difference in the outcrop size, where the PA site had a significantly larger outcrop area than the PB site.
3.1. Petrological Analyses and Physical Measurements
Characterization of the components for volcanoclastic deposits is usually limited to determining if they are monolithological or heterolithological, thus establishing this relationship if they can (or not) be related to an eruptive event [5] . Macroscopic analysis for both, volcanoclasts and the deposit’s matrix, were performed at the outcrops, and several clast samples were taken to measure their diameters.
Color was defined according to the Munsell Rock-Color Chart [28] , as well as texture, mineral components, alterations, and other physical characteristics to give the clasts a general compositional classification. The matrix grain size was obtained by directly measuring it with a sedimentation granulometer (sedimentometer), prepared from sediments coming from the Sacramento River in Chihuahua.
Physical properties were evaluated according to the following: 1) the mass of several samples was measured with an electronic weight scale; 2) the volume was obtained by immersing the samples in water inside a graduated container; and 3) the density of each sample was determined experimentally.
Since clasts have large sizes and the deposit’s height makes it difficult to sample, an optical technique via computer analysis was used to characterize the deposits. The method selected [29] had been already applied to volcanoclastic deposits. This method consists of two steps: 1) image analysis to select, quantify and measure clasts; and 2) functional stereology, which converts 2D data to a 3D environment using Perl code. This was done for one unit of the PA site deposit, while the other units and the PB site deposit were analyzed qualitatively as they did not contain the characteristics of area and number of clasts necessary for a quantitative analysis using the technique mentioned previously [29] .
Since it was observed how the quantity and size of the particles differed horizontally, a decision was made to work with two photos, taken at both ends of the unit. Although it is preferred to work with many samples as possible to discover trends among them [18] , the scale of the research did not allow such a detailed analysis.
3.2. Image Analysis
Image analysis is the part of the process where the clasts’ properties are quantified. It can be applied at several scales, from field photographs, or rock slabs’ scans, all the way down to thin sections. The collected image must be taken orthogonally to the sample or be digitally corrected. It must also contain an object of known length as a scale.
Comparing it to other methods, the sample’s scale does not influence the technique. However, it is necessary to consider how, both the area and the photograph’s resolution, limit the particle sizes to be measured: coarse particles will be surely affected by a small area, while a low resolution limits the visibility of finer particles [21] , although a low precision at both ends is inherent in the technique [29] .
Automated clast recognition by software is a fast and widely used technique in sediment analysis (e.g., [21] ), but it is difficult to apply in volcanoclastic deposits, since it rarely shows a distinctive contrast and boundaries between clasts [29] . In the case of the analyzed units, the in-between clasts contacts, a thin matrix layer in some of them, as well as a similar color range, caused too much “noise” for an automated procedure. So, clasts were outlined “manually” using the selection tools of an image editing program, which increased working time and possible errors. No distinction was made by clast type, given that they had a similar density. In addition, since it is an unconsolidated deposit, it was not necessary to consider diagenetic processes related to the clasts’ composition.
Once clasts were identified, images were analyzed using the ImageJ2 software included in the Fiji software package [30] , obtaining measurements of area, perimeter, particle shape descriptors (sphericity, aspect ratio, roundness/flatness), and Feret’s diameter. A threshold filter was applied to convert all color values to black and white, since false measurements when analyzing color images were displayed by the software.
3.3. Stereological Technique
Dimensions measured at a cross section are not an accurate representation of the actual measurements (such as clasts diameters) since the cut rarely coincides with the grain center [29] . Through the stereological techniques’ application, 3D knowledge of clasts can be obtained from the 2D observations of a random cross section [31] . Martin Jutzeler (University of Tasmania) and Alexander Proussevitch (University of New Hampshire), helped in this part of the process, by processing the 2D data obtained through their functional stereology technique. For details on the mathematical functions within the methodology, please consult the original source [29] .
The results of the stereological method are Gi values, which represent the “grain number distribution density”, equivalent to the number of particles per m2 or m3, and per class width (bin width), in units of (m2ϕ)−1 or (m3ϕ)−1, respectively. Separation between classes is ¼ϕ, since it offers greater precision when defining the distribution curve than the 1ϕ separation commonly used in sieving techniques [29] . From the Gi data, the volume fraction per class (Vi) is attained, which is the most common way of representing data obtained by optical techniques [4] .
However, to be able to compare with results coming from other studies (generally gained by sieving techniques), data must be converted to weight percentage (wt%) values. This can be done by multiplying Vi by the clasts’ density, whether this is an average, in case of having clasts of similar density (as in this case), or by the density of each type of clast.
The technique, although effective, truncates the analysis to coarser particles (<−1ϕ), since being performed on photographs, identification of smaller grains is practically impossible. It would be advisable to complement this information with a sieving analysis, however, this could not be done due to not having access to the sites. Therefore, the resulting data will have a certain “blind spot” in terms of the grain sizes identified.
For a statistical analysis of grain size data, it is routine to obtain certain parameters that represent a sum of differential thermogravimetric (DTG) characteristics [19] , such as: 1) average size (mean); 2) propagation (classification or sorting) of the sizes around the average; 3) symmetry or preferential propagation (asymmetry or skewness) towards one side of the average; and 4) degree of grain concentrations relative to the average (kurtosis), taken from the distribution histogram [32] .
It is typical in volcanoclastic deposit studies to obtain approximations to these parameters by employing graphical analysis, using the Inman formulations [4] [33] , or following other characterization methodologies [34] . On the other hand, in sedimentology studies there is a preference for the Folk & Ward techniques [18] [35] [36] , used to describe secondary volcanoclastic materials, while a combination of both techniques seems to be common for reworked deposits (e.g., [23] [37] [38] .
Following these techniques [23] , the GRADISTAT software [36] was used to perform the statistical analysis. This software uses the formulas in [35] and the “method of moments in statistics” to calculate the mentioned statistical parameters. The software code was modified to obtain the cumulative percentile values (D5, D16, D50, D84, D95) necessary to calculate the Inman parameters [33] .
4. Results
The sedimentological properties of mass flow deposits are important in distinguishing the phenomena [39] . Observation of relationships between the structures within the deposits located in both PA and PB sites (Figure 2), allows to obtain their general sedimentology parameters.
4.1. Field Description
PA site was exposed by an excavation pit of just over 6 meters deep and more than 50 meters wide. It was divided during field work into two blocks: the BNE and BSW Blocks (Figure 3). Each block has different characteristics, such as exposed surface and stratification, since BNE is well stratified, while BSW isn’t. A structure of 65˚/21˚ (azimuth and dip, respectively), serves a limit between both blocks.
![]()
Figure 2. Clast acquisition method applied to sample R1. Particles are smoothed during particle selection, as well as by corrections to avoid contact with each other.
![]()
Figure 3. PA site BNE and BSW blocks divided by an assumed fault, picture facing SE.
Although it was decided to define it as a fault, based on its own criteria, concrete evidence of the fault movement wasn’t found. Both blocks suffer an interruption in their upper part, being covered by a layer of unconsolidated material. The contact between this material and both volcaniclastic deposit blocks is flat and abrupt, suggesting an erosional process. This prevents analyzing the original volcanoclastic material topography, in addition to the intrinsic difficulties of a geological study in an urban area, which may be characteristic in its own definition [39] [40] .
Site B volcaniclastic deposits’ outcrop is small, with only 2 meters of visible material. Two types of facies can be distinguished: one with 4 to 30 cm clasts supported in a fine matrix, and another without clasts, showing a sand-sized granulometry (Figure 4). These units maintain a lower stratigraphic position compared to those of the PA site. They are also more consolidated, despite being exposed to weathering & erosion for a longer time. The larger clasts are quite fractured, although these fractures seem more a product of weathering than of transport.
Block BNE presents four visible units (or facies) (Figure 5), easily definable, 102˚/10˚ (azimuth and dip, respectively), each one with a different granulometry. From roof to base they are: I4) with a maximum thickness of 1.5 meters, composed of large portions of fine material forming a relatively massive layer with some gravel-sized clasts; I3) with a maximum thickness of 1.2 meters, composed of gravel-sized clasts and a fine light-colored matrix; I2) an intercalation of darker fine matrix with raw lamination and thin horizons of coarser grain size, barely exceeding sand size, and with the occasional presence of gravel-sized clasts; and I1) matrix-supported, with clasts similar in size to I3, but in smaller quantities and with some exceeding one meter in diameter (mega clasts).
4.2. Grain Analyses
In site PA, three clast compositions were recognized: 1) a tuff with intermediate composition lithics, 2) dacite and, 3) rhyolite. The first two seemed to be present in most, if not all, units within the BNE deposit, while only one rhyolite clast was
![]()
Figure 4. BSW units, aluvial soil layer covering the deposit.
![]()
Figure 5. BNE units de BNE. picture facing SE.
found within the transition zone (Figure 6). Adding “exotic” particles is common along the flow path, and most lahar deposits show a set of polylithological volcanoclasts [41] , a characteristic that also suggests a post-eruptive origin [5] .
According to Powers’ classification [42] , clasts show low sphericity and high roundness values (Table 1), possibly related to transport by water [43] , suggesting that most clasts are secondary [44] , something common in clay-poor lahars, originating as water flow which incorporates loose sediment from the volcano [41] . However, the exact percentage of clay in the deposit is unknown, and the high roundness values could well be due to the smoothing of the particles during image analysis (Figure 2).
![]()
Figure 6. Zone of horizontal change. It shows how I4 is maintained while I3 is interrupted, and similar grain-sized lenses appear. Rhyolite clast shown with a yellow line. View facing SE.
![]()
Table 1. Particle percentage within each roundness value, according to [43] .
An average density value of 2.63 g/cm3 was established for all clasts in general, based on measurements from 19 samples. This value lies in the lower part of the range reported for most igneous rocks, which supports the identification of the clasts as rocks of felsic to intermediate composition [45] and corresponds to the geology reported by the Mexican Geological Survey [27] .
Site PB clasts show a similar composition to those present in site PA. Although a 100% correspondence is unlikely, differences between both are only possible to be discerned through a microscopic study, so the presence of the same lithologies is assumed.
The finest material within the PA site was divided into two portions. The first was unconsolidated and supported the I3 clasts. It had an approximate thickness of 0.75 mm (very fine sand), grayish orange-pink color (Munsell 10R 8/2). The I4 unit appears to share this type of matrix, but this is based on mere observation, as it could not be analyzed. The second, found in units I1 and I2, is more compacted in certain sections. It has a moderately reddish orange color (Munsell 10R 6/6) and a coarser grain size of 0.149 mm (fine sand).
Two samples from unit I3 were analyzed. In total, 622 clasts from sample R1 and 181 clasts from sample R2 were measured, allowing comparisons between them and previous published work, since [18] [19] warn of the little usefulness of the quantitative study of a single sedimentary unit, and emphasize the comparison of DTG results from an unknown deposit with a known database.
The Gi output is equivalent to the number of particles per unit area or volume (depending on whether it is in 2D or 3D) and per size bandwidth phi. These results primarily serve as a first step to calculate volume per phi (Vi) and for a comparison between unit samples, since data from other deposits is rarely presented in units of (m2ϕ)−1 or (m3ϕ)−1. Both sample particles lie within the −1.5 and −7.75 ϕ (100% gravel) range, which although it may suggest a non-eruptive origin for the deposit, this may rather be due to the size truncation technique.
R1 has a unimodal distribution and, although a symmetrical distribution can be forced, there is no evidence that DTG behaves similarly. Meanwhile, R2 can be either monomodal or bimodal. Other work [41] reported lahars being commonly bimodal, especially in clast-supported facies, such as in I3. Also [39] , attribute bimodality to DAEVs and mention how it is eliminated when coarser particles deposit due to a debris flow, it being preserved in only some sections of the resulting deposit, this in the case of a flow originated by an AEV (air entry value). Although the unit’s bimodality is supported, it was decided to use the monomodal data from both samples to “standardize” the information (Figure 7).
Comparing the Vi distributions of the single-modal results of both R1 and R2, several differences can be inferred from both curves (Figure 8) and are described based on their statistical values (Table 2).
![]()
Table 2. Statistical parameters [33] [35] . Description from GRADISTAT [36] .
![]()
Figure 7. Grain size distribution (Gi) for R1 and R2 monomodals.
![]()
Figure 8. Volume size distribution comparison between R1 y R2.
First is how each sample has a different peak in its graph, with R2 mean value greater than R1. It was also noted how R2 covers a wider range of sizes (larger standard deviation) which is reflected in the slightly worse granular classification (sorting) compared to R1. This is opposed to lahars increasing their sorting with decreasing grain size, this, and the decrease in the asymmetry value with respect to the average grain size, coincides with data reported for DAEVs [39] . This increase in quantity and variety of particle sizes from R1 to R2 may indicate a possible flow direction [41] , in this case from SW to NE.
Mean size was plotted against σ1 to compare with different volcanic mass wasting flows (Figure 9). The closest position is of a non-cohesive debris flow (or clay poor), although it is necessary to consider whether this is not a result of the grain size distribution being cut at −1ϕ, resulting in much lower σ1 values.
When the samples are compared with data from lahars caused by rain, none is close to the R1 and R2 values; however, these show greater proximity to AEV values, specifically those of Río Pita [46] , and Cubilche [47] (Figure 10).
However, it should be noted that these values are plotted separately, without following the trend of the rest of their respective data, which becomes more visible when comparing the cumulative curves of both samples with those of Río Pita (Figure 11).
![]()
Figure 9. Comparison of R1 and R2 with [36] ranges for different types of mass wasting flows.
![]()
Figure 10. Comparison of R1 and R2 with data from lahars: Merapi [38] , and Semeru [48] , and AEVs: the Cubilche [46] , Cotopaxi [49] and Rio Pita [46] .
However, it should be noted that these values are plotted separately, without following the trend of the rest of their respective data, which becomes more visible when comparing the cumulative curves of both samples with those of Río Pita (Figure 11).
The definition of facies is made according to the type of deposition that occurs; therefore, it is necessary to establish this before continuing. Although the statistical evidence seems to indicate that the deposit corresponds to a volcanic debris avalanche (VDE), it lacks the weight to make a conclusion based on it. On the other hand, the rest of the physical evidence allows supporting the classification of the deposit as a lahar (Table 3).
![]()
Figure 11. Cumulative curve of R1 and R2 compared to information from Rio Pita [46] .
![]()
Table 3. Key characteristics of lahars and AEVs (modified from [39] [50] ).
The upper unit of the deposit (I4) has greater horizontal continuity, remaining without much change in the area where the rest is interrupted or mixed (Figure 5), except for the apparent incorporation of I3 clasts at its base. Therefore, it can be said that I4 belongs to a later phenomenon and little information can be obtained from its relationship with the rest of the units. Its fine granulometry, in the sand range (with the naked eye) and the isolated presence of clasts suggest that it is an over bank deposit [41] , but more information is needed to confirm this.
5. Conclusions
The genetic characterization of a volcaniclastic deposit located in a well know shopping mall in NW Chihuahua City, was performed by conducting field work, outcrop description as well as traditional sedimentological analysis as well as image analysis and functional stereology techniques.
The first site was classified as a secondary volcaniclastic deposit, product of two or more post-eruptive or non-eruptive lahars, defining two types of facies: 1) overbank facies (which overlies the rest of the units) and, 2) transitional facies with an inverse gradation, which disappears in a change zone which could be classified as lahar run out if a greater extent of the lahar were known.
A second outcrop was analyzed finding similarities in clast composition, but without any horizontal correspondence with the deposit. Even so, due to the proximity between these two, it can be assumed that it is the product of a similar deposition process.
It is recommended to expand the number of samples analyzed and to include petrographic analyses for both the clasts and matrix, as well as a granulometric analysis of the finest particles (<2 mm) to complement the grain size distribution and to obtain the clay percentage for I3 and the rest of the units.
In conclusion, this research has significantly contributed to the understanding of volcanoclastic sediments and has augmented the local geological knowledge base. The findings presented herein represent a foundational step for future endeavors, for academic pursuits, consulting engagements, or studies focused on identifying geological hazards within a city.
Acknowledgements
Foremost, our heartfelt gratitude extends to Drs. Martin Jutzeler of the University of Tasmania and Alexander Proussevitch of the University of New Hampshire for their invaluable guidance and unwavering support during the Stereological analysis. Their expertise significantly enriched the trajectory of our research. Furthermore, we express sincere appreciation to the Geology, Materials, Soils, and Asphalts Laboratories at the Facultad de Ingeniería (UACH) for generously providing the essential resources and steadfast support that proved instrumental in the successful completion of this project. Their collaborative efforts played a pivotal role in the realization of our research objectives.