Automated identification and mapping of geological folds in cross sections

2.4.1 Discrimination of composite folds and secondary folds

Composite folds mainly include two types: anticlinorium and synclinorium. A composite fold is a giant anticline (or a giant syncline) composed of different levels of secondary folds, which is the main structural style of an orogenic belt. They are also known as Alpino-type folds [31]. In the process of fold recognition, if a folded structure contains the core strata of other secondary folds, it can be identified as a composite fold. However, due to the nesting characteristics of composite folds, only the folded strata of secondary folds need to be processed in fold recognition and mapping. Therefore, in order to ensure efficient and reasonable fold model construction, composite and secondary folds within them should be correctly distinguished.

The general criterion of composite folds is as follows: if the stratigraphic interval of a fold overlaps with other fold intervals and contains its core strata, it can be judged as a composite fold; otherwise, it can be identified as a secondary fold (Figure 3). In addition, based on the interval scheduling method, the discrimination steps to distinguish between composite folds and secondary folds are as follows:

1. Read any two elements p _v , p _t (t ≠ v) from the set GD and judge whether the core position g _t of p _t is within the interval [e _v , r _v ]. If the result is true, then the element p _v is a composite fold, g _t is recorded as the attribute value of p _v .

2. According to equation (1), the position interval of p _v can be calculated as follows:

   (1) p v ∈ [ e v + r v − g t + 1 , g t − 1 ] if g v < g t , p v ∈ [ g t + 1 , e v + r v − g t − 1 ] if g v > g t ,

   where g _v is the core position of p _v .

3. If g _t is not within [e _v , r _v ], then the element p _v will be processed with other elements in GD, and the string s _v will be updated according to the final reduced position interval.

4. Traverse the set GD until the discrimination of the stratigraphic symbols of symmetrically repeated folds in GD is completed.

Figure 3

Schematic diagram for discriminating composite folds in cross section: (a) original cross section, (b) identification of the secondary folds in the cross section, (c) identification of the composite fold containing secondary folds in the cross section, and (d) extraction of the secondary folds at the core of the composite fold using the interval scheduling method.

2.4.2 Fold classification based on strike of the two limbs

According to the strike of the strata of the two limbs, the folded structure can be roughly divided into four types: isoclinal anticlinal fold, inclined anticlinal fold, isoclinal synclinal fold, and inclined synclinal fold [32] (Figure 4).

Figure 4

Schematic diagram of folds by a different strike in cross section: (a) isoclinal anticlinal fold, (b) inclined anticlinal fold, (c) isoclinal synclinal fold, and (d) inclined synclinal fold.

According to the stratigraphic strike information in the cross section, and through the new and old comparison between the core and its adjacent strata, the fold classification can be carried out. The specific steps are as follows:

1. Extract all stratigraphic line information in the section and store it in line set LN. Traverse all elements in the GD, query the corresponding stratigraphic line in LN, and get fold stratigraphic line set GL = {gl _v |v = 1,2,…,N} (Figure 5).

2. Read any folds p _v and find the corresponding stratigraphic line set gl _v . Read the stratigraphic line sl _v,f/2 and sl _v,f/2+1 and obtain the endpoints of the two stratigraphic lines: Q _a (x _a ,y _a ), Q _a′(x _a′, y _a′), and Q _b (x _b ,y _b ), Q _b′(x _b′, y _b′). Then distinguish the isoclinal fold and the inclined fold, according to the following equation:

   (2) Isoclinal fold if ( x a − x a ′ ) ( x b − x b ′ ) > 0 , Inclined fold if ( x a − x a ′ ) ( x b − x b ′ ) < 0 .

3. Read the central core character and adjacent character of the string s _v in p _v , and putting their stratum encoding into Astr and Bstr, respectively. The relationship between the new and old strata can be judged by comparing the value of their stratum encoding, which can be used to distinguish anticline and syncline (equation (3)):

   (3) Anticline if Astr > Bstr , Syncline if Astr < Bstr .

4. Traverse the set GD to judge all fold types based on strike.

Figure 5

Schematic diagram of a folded structure in cross section (line aa′, bb′, cc′, dd′, ee′, ff′: the stratigraphic lines; point a, b, c, d, e, f and point a′, b′, c′, d′, e′, f′: the endpoints of the stratigraphic line; fold stratigraphic line set: contains all stratigraphic lines within this folded structure).

2.4.3 Fold classification based on section morphology

Based on the section’s morphological characteristics, folded structures can be roughly divided into the following five types: gentle fold, open fold, close fold, tight fold, and isoclinal fold [20,33,34,35]. Its corresponding characteristics are shown in Table 2.

In cross sections, the shape of the hinge zone can be determined by calculating the interlimb angle and the fold angle (the angle between the normal vectors of two limbs) based on the attitudes of the stratigraphic line. According to equations (4) and (5), the fold angle α and interlimb angle ψ can be calculated, respectively. Furthermore, according to the classification standard in Table 2, the fold classification based on section morphology can be carried out (Figure 6).

where θ _a , θ _b are the dip angle of sl _v,a , sl _v,b . φ _a , φ _b are the inclination of sl _v,a , sl _v,b and sl _v,a , sl _v,b are any two corresponding stratigraphic lines in gl _v .

Figure 6

Geometry of a perfect fold in section view: (a) inclined fold , (b) isoclinal fold (i: the stratum-exposed location , h: the hinge zone , ψ: the interlimb angle , and α: the fold angle).

2.5.1 Boundary fitting method of two limbs of corresponding strata

After the classification of folded structures, the next step is to focus on fold mapping. Many studies have analyzed the shape of folds. For example, Stabler [36] and Hudleston [37] used Fourier analysis to analyze the shape of folds quantitatively; Bastida et al. [38] suggested fitting fold shapes to power functions. In this article, a simple scheme proposed by Srivastava and Lisle [39] and Srivastava et al. [40] was used to generate a smooth simulated folded curve based on Bézier curves. This method is easy to implement and modify, and the result is consistent with the natural fold shape. Furthermore, because the stratigraphic strike of the isoclinal fold and inclined fold differs, it is necessary to adopt different stratigraphic correlation mapping methods.

1. For an inclined synclinal fold (Table 3a(i)), the extension lines of stratigraphic lines in two limbs, which locate at point A and point B, respectively, meet at point C. Because the axial plane of the fold is also the angle bisector of the interlimb angle, the axial plane can be obtained by taking the angle bisector of ∠ACB, which intersects AB at point D. According to equations (6) and (7), the hinge zone endpoint E on CD can be found using Curve (the coefficient of the shape of the hinge zone):

   (6) Curve = DE CD ,

   (7) Curve = 0 – 0 . 3 if fold is gentle , 0 . 3 – 0 . 5 if fold is open , 0 . 5 – 0 . 7 if fold is close , 0.7 – 0 . 9 if fold is tight , 0.9 – 1 . 0 if fold is isoclinal ,

   where Curve is an interval range summarized according to the fold classification above, which can be used to express different fold section shapes. For different fold classifications, the parameter Curve needs to set different initial values (gentle fold: Curve = 0.2; open fold: Curve = 0.4; close fold: Curve = 0.6; tight fold: Curve = 0.8; isoclinal fold: Curve = 0.9).

   In order to ensure the smoothness and symmetry of the curve simulation of the two limbs of the fold, FG⊥CD through point E (Table 3a(iii)). Moreover, based on equation (8), Bézier curves are generated at points A, F, E and B, G, E.

   (8) x ( t ) = ( 1 − t ) 2 x a + 2 ( 1 − t ) t x f + t 2 x e , y ( t ) = ( 1 − t ) 2 y a + 2 ( 1 − t ) t y f + t 2 y e .

2. For an isoclinal fold, the axial plane is inclined, and the two limbs have the same tendency, so this fold needs to be discussed in two cases:

1. When the inclinations of the two limbs of the stratigraphic line are the same, but the dip is not the same, the correlation method between the two limbs is consistent with the calculation method of inclined folds (Table 3b).

2. When the two stratigraphic lines have the same dip and the same inclination, it means that the two stratigraphic lines are parallel, so we should calculate the formation thickness h between the stratigraphic lines and get the midpoint K between the surface endpoints A and B on the stratigraphic lines of the two limbs. According to the direction of the hinge zone, extend the distance h to point C from point K (Table 3c(ii)). Then get the endpoint E of the hinge zone by using Curve, where Curve = EK/CK. Finally, use quadratic Bézier curves to correlate the strata of the two limbs (Table 3c(iii)).

1. For a box fold or fan fold (Table 3d(i)), the dip of the two limbs is large, and the rock layer of the two limbs shows a partial overturn. The boundary fitting method of these fold types is similar to that of parallel stratigraphic lines (Table 3c) but the detail is differing. The midpoint K does not extend in the direction of the hinge zone but extends perpendicular to line AB. Moreover, the rest of the method is consistent.

Table 3

Association method based on Bézier curves

No.	Situation	Specific processing steps
		i	Ii	iii
a	Inclined fold
b	Isoclinal fold
c	Isoclinal fold with parallel stratigraphic lines
d	Box fold or fan fold

2.5.2 Optimizing the strata by the principle of superposition

In the process of single stratigraphic association, when the dip of the local stratum is too different from that of the adjacent stratum, the stratigraphic lines will cross (Figure 7a); when the folds are parallel folds, a single stratigraphic correlation cannot guarantee the consistency of stratigraphic thickness (Figure 7c). In order to solve these problems, an optimization method based on the principle of superposition was adopted.

Figure 7

The strategy of strata optimization treatment in cross section: (a) stratigraphic crossover due to a large difference in the dip, (b) optimized results of the anticline structure after adjusting Curve, (c) the formation thickness is not consistent in parallel folds, and (d) optimized results of parallel folds after adjusting Curve (point E _i : the endpoints of the hinge zones, point A′, B′, C′, D′, E′, F′: the endpoints of stratigraphic lines, h _i : the thicknesses of strata).

In nature, the adjacent strata have a certain similarity without being affected by tectonic activities such as faults. This similarity of adjacent strata is also called the principle of superposition [41]. The specific steps of the optimization method based on the principle of superposition are as follows:

1. Calculate the endpoints of the hinge zone of the folded curve. Traverse the smooth folded curves set GF, according to classification by folds, and calculate the endpoints E _i of all the folded curves.

2. Adjust the Curve according to the type of folds. When the fold is isoclinal, if the endpoint E _i of the outer stratum is lower than that of the inner stratum, the Curve of the outer stratum will be increased. Decrease the Curve of the inner layer successively to improve the endpoint of the hinge zone of the outer stratum (Figure 7b). When the fold is syncline, if the endpoint of the outer stratum is higher than that of the inner one, the processing procedure is the same as that of the anticline.

3. If the folds are parallel, the stratigraphic thickness of parallel folds should be consistent. Calculate the distance between all endpoints of the hinge zone E _i . According to the lithology, the distance between the hinge zone of the same strata should be kept consistent, which is consistent with the uniform thickness of the same strata in parallel folds (Figure 7d).

4. After optimizing all the folded curves, store them in the optimized fold curve set GF′.

2.6.1 Generating polygon features for folded strata

When the fold is synclinal, if the original stratigraphic line intersects within the section (Figure 8a), the original stratigraphic line corresponding to the synclinal fold will be replaced with the folded curves (Figure 8b).

Figure 8

Schematic diagram of fold formation in cross section: (a) when the original stratigraphic line intersects within the profile; (b) replacing stratigraphic lines with folded curves and generating polygon features of folded strata; (c) when the original stratigraphic line does not intersect within the profile; (d) combining stratigraphic lines with folded curves and generating folded strata (solid lines show stratigraphic lines and dashed lines show fold curves).

Suppose that the original stratigraphic line does not intersect the section (Figure 8c). In that case, the original stratigraphic line set LN will be combined with the folded curve set GF′ and construct the merged stratigraphic line set LF. The new stratigraphic lines and the surface lines of the cross-sectional construct planar strata. Furthermore, according to the original stratigraphic code and the stratigraphic code of the associated strata, attributes are assigned to the ground layer from left to right along the section direction (Figure 8d).

2.6.2 Optimizing folded strata by the principle of area conservation

A balanced cross section is one on which geometric criteria can restore structural deformation and displacement. Dahlstrom first explicitly proposed balanced cross sections in 1969 [42], and they follow three basic principles: volume conservation, area conservation, and layer-length conservation in a closed system. The theory of balanced cross sections is based on volume conservation; that is, the volume of a rock formation before and after deformation does not increase or decrease and remains constant. Under ideal conditions in a closed system, the stratum only changes its structure morphology during tectonic movement, but its volume remains unchanged [42,43]. However, in natural geological environments, it is impossible to reach the ideal conditions due to tectonic action, climate, deposition, denudation, and other factors in the evolution of the stratum. Therefore, the balanced cross-sectional technique is unsuitable for every geological section; it is an approximation based on basic geometric principles.

Considering that it is difficult to calculate the stratum volume in cross sections. This article adopts the theory of area conservation in balanced cross sections to optimize the strata. The specific steps of strata optimization are as follows:

1. According to the strata generated by the algorithm, calculate the area of any stratum S using Arc Engine API.

2. Calculate the length of the folded line and the stratigraphic line in this stratum and obtain the total length L (the length of this stratum before deformation).

3. According to the thickness of the same stratum, calculate the average thickness M of the stratum by using the following equation:

   (9) M = ∑ i = 0 n m i * 1 n ,

   where m _i are the layer thicknesses of the same stratum at various locations (e.g., m ₁, m ₂, m ₃, m ₄, m ₅, m ₆, m ₇ in Figure 9). m ₂, m ₄, m ₆ are the thickness of the hinge zone, which can be calculated by endpoints of the hinge zone E (see Section 2.5.1 for details). m ₁, m ₃, m ₅, m ₇ are the thickness of the stratum at the exposed position, which can be calculated by using the distance from the end point of the left stratigraphic line to the right stratigraphic line in cross section.

   (10) R = L * M .

4. Using the calculated total length L and average thickness M, calculate the pre-deformation stratum area R (equation (10));

5. Compare the stratum area S generated by the algorithm with the pre-deformation stratum area R. If S < R, it indicates that the area of the stratum model is too narrow. The narrower part of the stratum at the hinge zone (m ₆ in Figure 9) can be modified. The area of the stratum model can be increased by increasing the thickness of the stratum. If S > R, it indicates that the area of the stratum model is too broad. The wider part of the stratum at the hinge zone (m ₂ in Figure 9) can be modified. The area of the stratum model can be reduced by decreasing the thickness of the stratum.

Figure 9

Schematic diagram of a fold in cross section (L: the total length of the stratum; m _i : the layer thicknesses of the same stratum at various locations).

Of course, the area of the pre-deformation stratum obtained based on the above steps is only a rough calculation method and a simple estimation. A more reasonable optimization requires the support of more data.

The folded strata can be constructed and optimized through the above process, which effectively restores the folded structures in cross section and expresses some geological information. In the following sections, these processed cross sections are called “restored cross sections.”

3.2.1 Study area

The Parallel Fold Belt of Eastern Sichuan is located in southwest China, northeastern Sichuan Province, and the southern foot of Ta-pa Mountain, with a geographical range between 30°75′−32°07′ N and 106°94′−108°06′ E. Dazhou, with an area of 16,591 km² and an altitude of 222–2,458 m, is located in the Ridge and Valley Province of Eastern Sichuan. The general altitude trend is higher in the northeast and lower in the southwest. The geological structure of this area is diverse and complex. The southeast part belongs to the Eastern Sichuan Parallel Fold Belt of the Neocathaysian system, a typical isolated folded structure. The north and northeast parts belong to the arcuate fold belt in the southern margin of Ta-pa Mountain.

3.2.2 Experimental result

Draw three parallel section lines (Figure 11) along the vertical direction of the strata. The original sections without fold mapping were generated through the self-developed automatic section-drawing system (Figure 12a, c, and e). Using the method of this article, the folded structures in the original sections were identified and modeled successfully. Automatic restoration of the parsed cross sections (Figure 12b, d, and f) effectively shows the characteristics of the Parallel Fold Belt of Eastern Sichuan. It tells that, due to the lack of drill data, only the cross sections drawn by experts in the literature or geological map can be used to analyze the consistency of the final results. Since the corresponding cross sections cannot be found in this region, only a diagram of cross sections drawn by experts can be searched (Figure 13), so the section generated by the algorithm cannot be entirely consistent with it. However, the diagram of the Eastern Sichuan cross section tells that this area includes wide-spaced folds composed of a series of parallel closed narrow anticlines and gentle open synclines: anticlines become mountains, synclines become valleys, with valleys alternate and parallel to each other, which form the landscape of the parallel ridge valley area. By comparing with a diagram of a section drawn manually in this area, it can be found that cross sections created by the method of this article are consistent with the section diagram drawn by experts.

Figure 11

Parallel Fold Belt of Eastern Sichuan in the geological map.

Figure 12

Automatic drawing of the original Eastern Sichuan cross sections and the restored Eastern Sichuan cross sections: (a) an original Eastern Sichuan cross section (A–A′), (b) a restored Eastern Sichuan cross section (A–A′), (c) an original Eastern Sichuan cross section (B–B′), (d) a restored Eastern Sichuan cross section (B–B′), (e) an original Eastern Sichuan transverse cross section (C–C′), and (f) a restored Eastern Sichuan cross section (C–C′).

Figure 13

Diagram of the Eastern Sichuan cross section. (https://zhidao.baidu.com/question/2272689036589724068.html).

3.3.1 Study area

The Appalachian Mountains, extending along the eastern side of North America from Newfoundland (Canada) to Alabama (USA), are among the most recognizable and well-studied orogenic belts. The total length of the mountain chain is nearly 3,200 km, with a northeast-southwest trend. The southern Appalachian Mountains cover parts of five physio-graphic provinces. The Piedmont physio-graphic province is an upland of rolling hills with gentle slopes. The patterns of valleys and hills rarely coincide with the underlying bedrock structure. Most of the rocks of Piedmont formed as sediments or volcanic rocks on ocean floors, islands, and continental plates; igneous rocks formed when crustal plates collided, beginning about 450 million years ago. For this experiment, we selected the central Appalachian fold-thrust belt in the southern Appalachian Mountains. The central Appalachian fold-thrust belt has been described as a blind fold-thrust belt with few emergent faults [44,45], and the geographical range is between 37°15′−37°15′ N and 80°22′−80°30′ W (Figure 14).

Figure 14

Geological map of the southern Appalachian Mountains [46].

3.3.2 Experimental result

The original sections without fold mapping were generated through our self-developed automatic section-drawing system (Figure 15a and c). Using the method of this article, the folded structures in the original sections were identified and modeled (Figure 15b and d). Automatic restoration of the parsed cross sections effectively shows the characteristics of the central fold belt of the southern Appalachian Mountains. It tells that this area is a mesoscopic folding composed of a gentle and open anticline and a closed and narrow syncline [44]. Due to the lack of drill data, only the cross sections drawn by experts can be used to analyze the consistency of the final results. By comparing with a measured section (Figure 16) drawn manually on the geological map, it can be found that cross sections created by the method of this article are consistent with the measured section, and the overall stratigraphic strike and folded structure are consistent with the measured section.

Figure 15

Automatic drawing of the original Appalachian cross sections and the restored Appalachian cross sections: (a) an original Appalachian cross section (A–A′), (b) a restored Appalachian cross section (A–A′), (c) an original Appalachian cross section (B–B′), and (d) a restored Appalachian cross section (B–B′).

Figure 16

A cross section of the Appalachian Mountains drawn manually (C–C′) [46].

By analyzing the cross section of the Appalachian Mountains drawn manually (Figure 16), it can be found that the anticline is symmetrical above 300 m (the dip angle of the limbs is equal). However, below 300 m, the dip angle of the northern limb of the anticline suddenly steepens. The cross section automatically generated in this article is mainly drawn according to the attitudes of the outcrop, without considering the rapid change of strata strike below the surface. So, the modeling results are not inconsistent with the section drawn by experts. However, by comparison, it can also be found that the cross section generated by the method of this article is consistent with the measured real cross section in terms of overall stratigraphic strike and folded structure.

4.3.1 3D modeling of folded structures

A folded structure refers to the bending of rock strata under tectonic movement without losing its continuity and integrity. However, the strata of folded structures are incomplete due to denudation. This method of this article can not only construct the current denuded geological structures (Figure 17) but also support the reconstruction of denuded strata to restore the complete geological structures (Figure 18). After the reconstruction of the geological structure, the sections can be used to understand the spatial distribution of stratigraphic rock mass and to study and analyze the temporal and spatial variation of the stratigraphic rock mass. Based on the contrast between the 3D section model of the recovered folded structure (Figure 18) and the model of the incomplete folded structure (Figure 17), we can infer the temporal and spatial evolution process of this area or analyze the influence of the crustal movement.

Figure 17

3D model of unrestored structures constructed from Eastern Sichuan cross sections.

Figure 18

3D model of the recovered folded structure constructed from Eastern Sichuan cross sections.

4.3.2 Mapping of folded structures

In addition, different from the traditional cross section, this algorithm can construct folded structural elements in the section (Figure 19) and express the attribute information (e.g., fold type, fold morphology, fold scale) and specific ranges of geological structures (Table 4).

Figure 19

Mapping and expression of folded structural elements in the Eastern Sichuan cross section (C-C′) (AF: anticlinal fold; SF: synclinal fold).

This article mainly infers and correlates the folded strata below the surface in the cross section according to the stratigraphic morphological characteristics of folded structures. Since the reconstruction of the eroded part of anticline structures is not the focus of this article, only a rough method is provided here. In order to reconstruct the eroded part of the anticline accurately, the mechanism of folding and the structure state of the layers must be considered. For example, folded competent layers should be thinned toward the crest of the anticline if they are as growth strata (as syn-tectonic layers). If the same layers are pre-tectonic, they should not change in thickness toward the crest of the anticline. For restoration and reconstruction of layers before deformation, examining the folding mechanism and how it has deformed the layers is necessary. In the future, more detailed inference rules can be considered to more precisely model the eroded part of the anticline.