1. Introduction
Modelling of geological objects is a prerequisite for reconstructing the evolution of the geometry of a structure through time. Geological structures are by essence known from sparse and heterogeneous data, excluding 3D seismic cubes, and are generally difficult to interpret owing to the lack of data and to the inherent complexity of these objects. Moreover, the available data are heterogeneous: a model has to take into account intersections of lines and planes (wells and horizons for example), both 1D, 2D and 3D data from maps, wells and seismic acquisitions, vectors and tensors (velocity, stress and strain fields) and rock physical properties. These data are also irregularly distributed and generally clustered, and show various degrees of confidence (well and field data are generally ‘harder’ than seismic data) and have different needs for interpretation and processing (time to depth conversion of seismic data). Numerical methods known as computer-aided design have been widely used for the modelling of natural or manufactured objects, such as car designs, based on continuous mathematical description of physical space (Bézier-like curves, see Bézier, Reference Bézier, Barnhill and Riesenfeld1974; Farin, Reference Farin and Farin1988). However, the major drawback of these early methods is their continuity regarding the discrete character of the geological data. Although several authors proposed approaches to overcome this problem, a major breakthrough was made in the development of discrete methods based on a meshing of the surfaces and volumes similarly to the finite element techniques (see Mallet, Reference Mallet2001 and numerous references therein, and Caumon, Reference Caumon2009). Based on a set of mathematical methods allowing modelling of the topology, the geometry and the physical properties in a discrete space, geomodelling, as defined by Mallet (Reference Mallet2001) and now widely used in numerical tools such as Gocad or Petrel, is a powerful method to describe geological space in a satisfactory manner (e.g. Caumon, Reference Caumon2009).
Geomodelling of geological structures has many applications (3D localization, volumetric estimates), and in particular it provides a starting point for the process of geological structure restoration, an operation through which the initial state (i.e. pre-deformational state) of a geological object can be built. In the restored state, one of the most common assumptions is then that the different horizons should have equal properties, such as length or area, so-called balanced (Chamberlain, Reference Chamberlain1910). This technique has been used by numerous authors since (Buchner, Reference Buchner1933; Goguel, Reference Goguel1962; Laubscher, Reference Laubscher1962). Balanced sections have been generalized since 1969 (Dahlstrom, Reference Dahlstrom1969; Suppe, Reference Suppe1983), and progressively incorporated into numerical tools (Moretti & Larrère, Reference Moretti and Larrère1989; Rowan & Kliegfield, Reference Rowan and Kligfield1989). Following the development of acquisition techniques from 2D to 3D, the attention of structural geologists focused on 3D data. Following the development of 2D cross-section balancing tools, a step further was reached by developing map restoration (Gratier et al. Reference Gratier, Guillier, Delorme and Odonne1991; B. Guillier, unpub. Ph.D. thesis, Univ. Joseph Fourier, Grenoble, 1991; Rouby et al. Reference Rouby, Cobbold, Szatmari, Demercian, Coelho and Rici1993; Rouby, Xiao & Suppe, Reference Rouby, Xiao and Suppe2000; Dunbar & Cook, Reference Dunbar and Cook2003). The techniques used to restore maps range from pure geometrical methods (block mosaic, e.g. Gratier et al. Reference Gratier, Guillier, Delorme and Odonne1991) to solutions of mechanical problems where horizons are treated as hyperelastic membranes (Dunbar & Cook, Reference Dunbar and Cook2003). Nevertheless, even if the maps are treated as data within the 3D space, these types of restoration are not really 3D. Real 3D restoration that handles 3D data as such became only recently available (P. Muron, unpub. Ph.D. thesis, Institut National Polytechnique de Lorraine, Nancy, 2006; Moretti, Lepage & Guiton, Reference Moretti, Lepage and Guiton2006; Maerten & Maerten, Reference Maerten and Maerten2006; Durand-Riard, Caumon & Muron, Reference Durand-Riard, Caumon and Muron2010).
Following the principles of 3D geological modelling, the objectives of this structural study of the Zagros fold-and-thrust belt were to assimilate and exploit a vast set of data acquired by successive generations of geoscientists ranging in age from the 1960s to the present day. This included geological maps, field observations of structures, reports, logs and checkshots of exploration wells, reports and cross-sections drawn for individual structures, and 2D seismic lines. To this end, a 3D digital geological model has been constructed using a geomodeller, both to express ideas concerning the regional structures and to integrate a large amount of data from all the sources previously mentioned in a coherent manner. The Lurestan area is located in the southern part of the Zagros fold-and-thrust belt, in the southern Fars province of Iran, and constitutes the target of this study. In this part of the belt, also known as the Fars arc, the presence of multiple pre-existing salt diapirs developed since Early Palaeozoic time complicates the fold geometry and explains the great difficulty in extrapolating geometry at depth when subsurface data are lacking (Callot, Jahani & Letouzey, Reference Callot, Jahani, Letouzey, Lacombe, Lavé, Verges and Roure2007; Jahani et al. Reference Jahani, Callot, Frizon de Lamotte, Letouzey, Leturmy, Lacombe, Lavé, Verges and Roure2007, Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009; Callot et al. in press). The studied area and the structural model cover around 60000 km2, taking into account over 70 anticline structures and 73 emergent Hormuz salt diapirs (Fig. 1a). The construction of the structural model confronted many of the technical and geological problems to be faced in the southern Fars area, namely taking into account over a considerable surface area numerous and closely spaced fold structures, with coexisting emergent and buried salt diapirs, variable thicknesses and sedimentary facies of the strata, sparse seismic data of generally low quality, and sparse well control.
(a) Geological map of the southern Fars (Bandar Abbas area), Iran (from NIOC, 1:100000, 1969). a – Hormuz salt diapirs, Infracambrian to Lower Cambrian; b – Upper Triassic–Jurassic–Lower Cretaceous; c – Middle Cretaceous; d – Upper Cretaceous; e – Paleocene–Eocene–Oligocene; f – Lower–Middle Miocene; g – Upper Miocene; h – Late Upper Miocene–Lower Pliocene; i – Upper Pliocene–Quaternary (modified from Jahani et al. Reference Jahani, Callot, Frizon de Lamotte, Letouzey, Leturmy, Lacombe, Lavé, Verges and Roure2007). The outline of the study area is shown in pink (grey). Squares refer to the following figures: A –Fig. 10; B – Fig. 8; C – Fig. 12; D – Fig. 11. (b) Simplified lithostratigraphic chart of the eastern Zagros fold-and-thrust belt and southeastern Persian Gulf.
2. Geological setting
2.a. Tectonostratigraphic setting
The Zagros fold-and-thrust belt belongs to the Alpine–Himalayan orogenic belt, and extends NW–SE over 1800 km from the Taurus Mountains (northeastern Turkey) to the Hormuz Strait, Iran (Stocklin Reference Stocklin, Burk and Drake1974; Haynes & McQuilan Reference Haynes and McQuillan1974; Alavi, Reference Alavi1994). The Main Zagros Thrust (MZT) constitutes its northern limit and it dies out in the Persian Gulf, which represents its present-day active foreland basin. The eastern part of the Fars province constitutes the studied zone (Fig. 1a).
The following stages of structural evolution are generally recognized in the Zagros fold belt: (1) a Palaeozoic platform phase; (2) a Permo-Triassic rifting phase which only affected the northern and eastern boundaries of the Zagros platform; (3) Jurassic to Early Cretaceous platform and passive continental margin of the Neo-Tethys Ocean; (4) Late Cretaceous ophiolite emplacement (obduction), and finally, (5) following the beginning of the collision phase dated back to Eocene time (Fakhari et al. Reference Fakhari, Axen, Horton, Hassanzadeh and Amini2008; Ballato et al. Reference Ballato, Uba, Landgraf, Strecker, Sudo, Stockli, Friedrich and Tabatabaei2011; Mouthereau, Reference Mouthereau2011, this issue), in accordance with detrital zircon ages of 35 Ma (Allen & Armstrong, Reference Allen and Armstrong2008; Horton et al. Reference Horton, Hassanzadeh, Stockli, Axen, Gillis, Guest, Amini, Fakhari, Zamanzadeh and Grove2008), a Neogene phase of exhumation, folding and crustal shortening (e.g. Agard et al. Reference Agard, Omrani, Jolivet and Mouthereau2005; Mouthereau et al. Reference Mouthereau, Tensi, Bellahsen, Lacombe, Boisgrollier and Kargar2007). Salt tectonics probably began in Early Palaeozoic time (Jahani et al. Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009), and actual diapir activity has occurred at least since Permian time, and has continued up to the present day (R. A. Player, unpub. Ph.D. thesis, Univ. Reading, 1969; Ala, Reference Ala1974; Mottei, Reference Motiei1995; Jahani et al. Reference Jahani, Callot, Frizon de Lamotte, Letouzey, Leturmy, Lacombe, Lavé, Verges and Roure2007, Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009).
During Early Palaeozoic time (Fig. 1b), shallow marine to fluvial clastic sediments extended onto the erosion surface formed either on the Precambrian basement or on the pre-Hormuz basins. The whole area is transgressed by shallow marine sediments deposited on early clastic sediments, following the regional Hercynian unconformity (Faraghun Formation; Berberian & King, Reference Berberian and King1981). The Late Permian and Early Triassic periods in the Zagros correspond to an episode of marine carbonate sedimentation followed by an evaporitic platform (Murris, Reference Murris1980; Koop & Stoneley, Reference Koop and Stoneley1982; Setudehnia, Reference Setudehnia1978), laterally evolving to dolomite to the north (Setudehnia, Reference Setudehnia1978; Szabo & Kheradpir, Reference Szabo and Kheradpir1978; Murris, Reference Murris1980).
During Jurassic to Middle Cretaceous times (Fig. 1b), sediments were deposited in a basin controlled by vertical movements and flexures along major basement faults (Berberian & King, Reference Berberian and King1981). Shortening of the Arabian margin began on its northeastern edge during Early Coniacian–Late Santonian times (Ricou, Reference Ricou1971; Falcon Reference Falcon and Spencer1974). Foreland deposits of Campanian to Maastrichtian age (deep water marls, shales, marly limestones and turbidites) were deposited in front of the thrusted ophiolite and basinal/slope pile (Homke et al. Reference Homke, Vergés, Garcés, Emami and Karpuz2004, Reference Homke, Verges, Serra-Kiel, Bernaola, Sharp, Garces, Monteru-Verdo, Karpuz and Goodarzi2009). The Upper Cretaceous and Palaeogene sediments display highly variable thicknesses and facies owing to their progressive deformation during the ophiolite obduction and early collision (Koop & Stoneley, Reference Koop and Stoneley1982; Alavi, Reference Alavi1994; Agard et al. Reference Agard, Omrani, Jolivet and Mouthereau2005; Horton et al. Reference Horton, Hassanzadeh, Stockli, Axen, Gillis, Guest, Amini, Fakhari, Zamanzadeh and Grove2008; Homke et al. Reference Homke, Verges, Serra-Kiel, Bernaola, Sharp, Garces, Monteru-Verdo, Karpuz and Goodarzi2009).
Continental collision started in Late Eocene time and propagated during Oligocene time at the northern promontory of the Arabian plate (Agard et al. Reference Agard, Omrani, Jolivet and Mouthereau2005; Horton et al. Reference Horton, Hassanzadeh, Stockli, Axen, Gillis, Guest, Amini, Fakhari, Zamanzadeh and Grove2008), and propagated southeastwards during Early Miocene time (early growth strata in Agha Jari sediments aged 15 Ma; Khadivi et al. Reference Khadivi, Mouthereau, Larrasoaña, Vergés, Lacombe, Khademi, Beamud, Melinte-Dobrinescu and Suc2010). The main phase of folding is recorded in the Upper Agha Jari Formation (Upper Miocene–Pliocene) by growth strata (Berberian & King, Reference Berberian and King1981; Homke et al. Reference Homke, Vergés, Garcés, Emami and Karpuz2004; Sherkati, Letouzey & Frizon de Lamotte, Reference Sherkati, Letouzey and Frizon de Lamotte2006; Khadivi et al. Reference Khadivi, Mouthereau, Larrasoaña, Vergés, Lacombe, Khademi, Beamud, Melinte-Dobrinescu and Suc2010). Nevertheless, there is evidence of earlier tectonic movements recorded in the Asmari and Gachsaran formations (Hessami, Koyi & Talbot, Reference Hessami, Koyi and Talbot2001; Sherkati, Letouzey & Frizon de Lamotte, Reference Sherkati, Letouzey and Frizon de Lamotte2006). Eventually, the Bakhtyari Conglomerate (Plio-Pleistocene) followed the main phase of folding, although it is locally affected by ongoing deformation. Recent kinematic scenarios proposed for the Zagros fold-and-thrust belt (Molinaro et al. Reference Molinaro, Leturmy, Guezou, Frizon de Lamotte and Eshragi2005; Sherkati, Letouzey & Frizon de Lamotte, Reference Sherkati, Letouzey and Frizon de Lamotte2006) suggest a two-step evolution with a mainly thin-skinned phase during Miocene time followed by a basement-involved phase since Pliocene time (Fig. 1b). However, the fact that both cover and basement are currently deforming (Roustaei et al. Reference Roustaei, Nissen, Abassi, Gholamzadeh, Ghorashi, Tatar, Yamini- Fard, Bergman, Jackson and Parsons2010), and that several lines of evidence support the early involvement of the basement in shortening during the Zagros history (Lacombe et al. Reference Lacombe, Mouthereau, Kargar and Meyer2006; Mouthereau, Lacombe & Meyer, Reference Mouthereau, Lacombe and Meyer2006; Mouthereau et al. Reference Mouthereau, Tensi, Bellahsen, Lacombe, Boisgrollier and Kargar2007; Ahmadhadi, Lacombe & Daniel, Reference Ahmadhadi, Lacombe, Daniel, Lacombe, Lavé, Verges and Roure2007), indicates that in the Fars province, both cover and basement were instead more or less coevally shortened through superimposed thin-skinned and thick-skinned tectonic styles of deformation.
2.b. Halokinesis
The Hormuz salt layer pinches out southwestward against the normal palaeofaults bounding the Qatar Arch and the stable Arabian plate, but it extends northward in the Persian Gulf and southeastern Fars province. In the eastern Zagros it has been proposed that salt movement largely preceded folding (Harrison, Reference Harrison1930; Kent, Reference Kent1958, Reference Kent1970; R. A. Player, unpub. Ph.D. thesis, Univ. Reading, 1969; Jahani et al. Reference Jahani, Callot, Frizon de Lamotte, Letouzey, Leturmy, Lacombe, Lavé, Verges and Roure2007). These propositions rest on field observations, namely of reworked Hormuz material and progressive unconformities (i.e. growth strata) in the vicinity of salt diapirs, and with seismic evidence of downbuilding of lateral mini-basins or rim synclines around salt domes (Jackson & Talbot, Reference Jackson and Talbot1991; Rowan, Jackson & Trudgill, Reference Rowan, Jackson and Trudgill1999; Rowan & Vendeville, Reference Rowan and Vendeville2006; Al-Barwani & McKlay, Reference Al-Barwani and McKlay2008; Jahani et al. Reference Jahani, Callot, Frizon de Lamotte, Letouzey, Leturmy, Lacombe, Lavé, Verges and Roure2007, Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009). The salt movement has occurred in the southern Fars since Late Cretaceous time for many diapirs, and for some of them even since Late Permian time (Motiei, Reference Motiei1995). But from seismic evidence, salt doming and dome growth started in Early Palaeozoic time and continued up to the present day (Jahani et al. Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009).
During Early Palaeozoic times mini-basins were formed as a result of clastic sediment progradation with syn-sedimentary salt movement resulting from differential loading. Pre-existing basement structures most probably controlled the orientation and location of the mini-basins above the salt and associated salt ridges, which evolved quickly towards nearly circular salt bodies sometimes surrounded by weld during salt withdrawal, resulting in buried diapirs. On the contrary, diapirs still fed by the salt layer continued to grow during Late Palaeozoic–Mesozoic times, and stayed close to sea level with episodic emergence and erosion. The Tertiary compressive events reactivated the ascent of salt for most of the domes, and a majority of them breached the surface owing to erosional processes during folding. This resulted in a dramatic increase in salt diapir diameter above the Asmari Formation in a characteristic ‘mushroom’ shape (Callot et al. 2011).
3. Data
3.a. Dataset
The National Iran Oil Company (NIOC) compiled for this study a vast dataset in digital format (see Fig. 2 for locations of the available seismic lines and wells). The NIOC made a considerable effort to gather the data from various sources, and in a variety of formats. In particular, a set of 36 well data was received from both inside and in the close vicinity of the project area. The level of data and detail varies from well to well, according to the history of exploration in that location, and the preserved data.
Basemap of the study area wrapped on topography, illustrating the seismic (in paper or electronic format) and well data base available for the study.
Well summary sheets were available for all wells, describing layers encountered, stratigraphic marker depths and well results. For some of the wells, log data in LAS (Log ASCII Standard) and image formats were provided, that allowed us to check the pertinence of the summary sheets. Palaeologs of biostratigraphic analyses were seldom available. Finally, some checkshot measurements were supplied for a few wells.
In addition to well data, information was recovered on 70 individual structures throughout the study area and close to its limits, mainly as reports and/or extracts of reports in English and Persian languages, but also consisting of structural cross-sections made during different periods and by different geologists with their location map. Over a hundred of these sections were incorporated in the model.
Outcrop stratigraphic logs were also gathered. The database was also supported by field and helicopter surveys on the structures. The surveys included the vast majority of the structures and the diapirs of the area. The observations were taken into account in the geomodel concerning the general structural style, the coherency of dips between the various layers, the disposition of salt plug outlines, and to control and verify the emergence of fault planes at the surface when the faults were modelled. Most wells are drilled through Mesozoic strata and only a few encountered Upper Palaeozoic rocks.
The precise locations of the structures and the identification of the layers cropping out were supported by recent 1:100000 scale geological maps of good quality, that were available for the whole of the study area with the exception of one missing quadrant sheet (see Section 4.d). Dip and strike data for the Top Asmari surface was used to control the geomodel. The intersection of the Top Asmari surface with topography was also used as a control.
Finally, approximately 130 seismic lines were available for interpretation and incorporation into the model. The seismic data was not homogeneous, consisting of SEG-Y data and scanned seismic data (Fig. 2).
3.b. Data classification
Because of this variety of sources, for modelling purposes we separated the information into ‘hard’ data and ‘soft’ or interpreted data. The hard data should be considered an immutable basis for the modelling, and all new models applied to the area must respect these data. On the contrary, soft or interpreted data are subject to possible changes, and any new interpretation will necessarily impact the image of the model.
Hard data for the modelling resulted from direct measurements or values that are undisputed and approved. We can classify in that category in particular, surface data, as the topographic DEM values, but also surface geology observations, either obtained from geological maps or from direct field observations of the structures. The stratigraphic layer identification or outcrop stratigraphic logs, and also the dip and strike values of the layers were included in that category. We also classify in the hard data category the stratigraphic marker depths obtained in the wells of the area, as this data was considered highly reliable, supported by well log data and biostratigraphic analysis from well samples. By deduction the checkshot values, which allow the deduction of the formation velocities for some wells of the area, were also considered hard data at the well locations.
It should be noted that the surface geology mainly gives information on the Tertiary and Mesozoic layers from Top to Base Cretaceous age (Fig. 1a), and with very few older outcrops. The only Palaeozoic formations cropping out that cover ages down to Ordovician are restricted only to the extreme northeast of the study area, and solely on the structural trend of the Kuh-e-Gahkum and Kuh-e-Faraghan area. The control of the deep layers is therefore not well provided by outcrops.
3.c. Seismic information
The main challenge of the 2D interpretation was the patchy nature of the data coverage: small seismic campaigns shot mostly in the 1970s had been used historically for spot evaluation of single structures or groups of structures. The data is located mostly in the flat, intermontane plain areas of the study area. The different seismic campaigns were at datum planes ranging from sea level to 1000 m above sea level. Additional small shifts had to be applied to allow tie-in of the different campaigns at crossing points.
Soft/interpreted data mainly consist of seismic sections and structural underground cross-sections covering the study area. Although seismic information consists of physical measurements it should be regarded as soft data for several reasons. First, the seismic data output is in the time domain, and the model deals with depth representation. There is therefore some uncertainty attached to the conversion from time to depth.
Secondly, a large number of the seismic lines in the data package were shot prior to the 1980s. The quality is moderate to poor. Good sections are related to synclines, or in areas with low dips of the layers. However, sections covering high dips are generally of poorer quality. In addition, there are static problems, especially for structures rising well above mean sea level. For these structures, the shallow horizons are often impossible to pick with high confidence.
Thirdly, the sections available were processed in the 1970s, and generally stacked and not migrated. Reprocessing of a few lines where raw data exists may improve locally the precision of the model.
Finally, the different surveys are frequently spatially disconnected from each other, with no direct ties (Fig. 2). The identification of non-cropping out horizons therefore relied on well control where available, but this was practically never the case, and hardly ever with a synthetic seismogram. In other cases, tying between lines relied on regional modelling of time isopachs and direct correlation of seismic character.
For all these reasons, it must be considered that the identification of the seismic horizons mapped is subject to reinterpretation rather than being a hard ‘input’ of the modelling process.
4. Methods
4.a. Introduction
The construction of 3D geological models from sparse data has been explored in detail by geomodelling software development teams over a period of several years. In particular, specific techniques have been developed to construct a surface which honours data points using interactive editing of triangulated surfaces (e.g. the Gocad software, see Caumon, Antoine & Tertois, Reference Caumon, Antoine and Tertois2007; Franck, Tertois & Mallet, Reference Franck, Tertois and Mallet2007).
We chose to use the wealth of information brought by the field data contained on the geological maps and from observations made during the field surveys to supplement the much sparser control provided by seismic lines and wells.
The modelling strategy (Fig. 3) was to first load into a geomodeller and in the depth domain, both topographic and geological data (DTM, 1:100000 scale geological maps, satellite images, wells and stratigraphic marker depths). Then, interval thickness data from wells both inside and outside the study area were incorporated into isopach maps and used for velocity analysis to prepare the integration of depth-converted seismic picking.
Modelling strategy: A detailed geological model was built of one key surface, which was then used in turn to construct all of the underlying surfaces using a regional depth isopach model fed from both seismic and well control.
The seismic data obviously had to be converted from time to depth before being used in the geomodel. Numerous conflicts between the near-surface data already entered into the geomodel and the seismic pick required the iterative refinement and re-picking of the seismic horizon and adjustment of the regional velocity model in order to honour surface data and well picks. At the end, a fairly simple but robust velocity model was used to convert the seismic horizon to depth.
4.b. Seismic velocity conversion
Velocity data at wells (checkshots) were used to create a regional velocity model, but the seismic velocity conversion model is still debated. S. Jahani (unpub. Ph.D. thesis, Univ. de Cergy Pontoise, 2008) already addressed in detail the question of a regional velocity model for the southern Fars and offshore Persian Gulf regions, with far more wells than were available in the scope of the present study. However, some modifications had to be applied to this regional model in order to better fit seismic data in the onshore area: a simple model of one interval velocity per layer based uniquely on the onshore well checkshot data proved to perform better than either a velocity model derived from onshore and offshore wells or (surprisingly) than velocity maps interpolated between well control points (Fig. 4).
Comparison of the adequate fit of the seismic picking with horizons deduced from hard data constraints and converted to time with various time-depth models. Solid lines are the result of seismic picking. Dashed lines are result of time-converted geomodel horizons. When solid and dashed lines superpose, the model is better. Left-hand side: constant velocity per layer model integrating data from onshore wells and offshore wells. Central panel: constant velocity per layer model integrating data from onshore wells only (retained model). Right-hand side: variable velocity per layer model, obtained by interpolation of the well data.
Firstly, the results of interval velocity calculation for the onshore wells based on checkshots were compared to the results of a more regional survey of well interval velocity data. It appeared that the offshore wells included in the regional survey have the effect of decreasing the calculated regional average velocities. This is hardly surprising given the presumed greater palaeo-burial depth and compaction of the series onshore compared with the relatively undisturbed offshore sediments. An average based only on onshore wells is therefore logically more representative of our onshore study area.
In addition, the regional velocity maps presented by S. Jahani (unpub. Ph.D. thesis, Univ. de Cergy Pontoise, 2008) and re-created for the sensitivity tests contained some control points with extreme values (for example, for the Guniz-1 well, velocities for the shallow level are particularly high). Our interpretation is that they probably correspond to fairly local velocity anomalies, although the cause of such high values is still to be investigated. It seemed, however, doubtful that these latter values should be used to create a regional velocity model.
As a general rule, it was felt that interpolation between such sparse data points results in a velocity map which is largely under-constrained outside of the control points compared to a simple linear velocity model. As a result, the velocity model was kept constant for each interval and the isopach model was changed in order to create a best-fit set of data in depth and time. To account for precise fit to wells, it should be noted that the final step of the modelling of the horizons included a last local fit to well marker data, which can ultimately be considered to be a form of local correction of the inaccuracy of the simplistic regional velocity model.
4.c. Maximizing hard data integration
In order to optimize the integration of hard data in the model, we decided to build the detailed geological model based on one key surface, for which controls could be maximized. This surface is in turn used as a reference surface to construct underlying surfaces using a regional depth isopach model fed from both seismic and well control. The Top Oligo-Miocene Asmari horizon was a natural choice for this key horizon, for several reasons. First, the Asmari carbonate presents a very characteristic facies, which in most of the study area marks the limit between the pre-Miocene collision folding and the post-Miocene kinematics. It is a shallow horizon, cropping out over most of the study area, therefore providing frequent control points where its intersection with topography and dip/strike data can be integrated into the model.
Secondly, our field observations showed that dips measured in horizons above the Asmari Fm were often decoupled from the structure below, because some of them are syn-tectonic (growth strata for example), but also through the influence of a mobile intermediate horizon immediately situated above the Asmari Fm (the Gachsaran evaporite and lateral equivalent), making it impossible to model deeper horizons from dip information taken from horizons shallower than the Asmari Fm.
Thirdly, the Top Asmari is a regional seismic horizon, being both relatively easy to follow along and between seismic sections, and frequently tied in with outcrop data. The tie is generally well documented, because the Asmari Fm is systematically sampled by wells as a regional petroleum reservoir, and fairly reliably identified. However, the Top Asmari cannot be considered to be a reliable time line. The Asmari carbonates interfinger with the Razak Formation in the northeast of the study area, and are locally undifferentiated from the underlying Jahrum carbonates (Jahani et al. Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009). The Top Asmari was therefore modelled as a major lithological boundary, not as a time surface.
Once a reliable model for the Top Asmari surface had been built in depth, regional isopach maps for the key intervals below were used to construct deeper and deeper surfaces stepwise down to a phantom deep horizon taken as 2000 m below the Faraghun seismic marker. Obviously this later Lower Palaeozoic interval merely allows the creation of a minimum estimate for the real depth of the Hormuz salt horizon.
4.d. Controls on the Top Asmari surface model
As mentioned previously, the available hard data which can be used to control the geological model of the Asmari surface are the traces of the outcrop of the Asmari limestone horizon, field-measured dips and strikes and well marker information. The contour of the cropping out Top Asmari surface was traced from a mosaic of 1:100000 scale geological maps (unpublished geological maps from the NIOC exploration directorate: 20857E, 20858W, 20858E, 20859W 20859E, 20862E, 20863W, 20863E, 20864W, 20864E, Orzuieh, 20868E, 20869W, 20869E 20870W, 20870E, 20871W, 20874E, 20875W, 20875E, 20876W, 20876E, Chiru, Charak, Band-e Lengeh, Dojgan, 20881E, according to internal classification of the NIOC EXP) using satellite images as a control, which was especially important in the rare cases where the mapped contours failed to meet at the boundary between different maps. At the same time, control points for bedding dip and strike were digitized wherever available on the maps. Well markers were taken directly from existing databases and fed into the model (Fig. 5a).
(a) First modelling step for Asmari surface modelling, using only hard data (outcrop of Asmari, field dip and strike directions and Asmari well markers). (b) Second modelling step for Asmari surface modelling, with integration of soft data (depth-converted seismic picks for the synclines and interpreted geological cross-sections for structures above erosion level). The resulting model of the Asmari surface had a grid size of around 5 km in a structure-parallel direction and 300 m in a structure-perpendicular direction. Residual errors in well marker fit ranged from 0 to 25 m, essentially owing to the difficulty of honouring a closely spaced control point array with a still relatively coarse mesh.
The soft data, available in the form of depth-converted seismic picks, was incorporated to constrain each horizon's maximal depth in the synclines, as we cannot extrapolate precisely the depths in the synclines from surface data. As said previously, most of the recent formations are syn-kinematic, and observed dips and thicknesses at the surface may lead to underestimation of depths in the core of the synclines.
We used also geological cross-sections to model the Top Asmari above and near the topographic surface because they represent a geometrical construction or extrapolation from field observations along the section line. The hundreds of detailed geological cross-sections, which reflected the work carried out by various geologists in numerous field campaigns, contained most important detailed observations of bedding dips and strikes and local formation thicknesses. These data were a useful addition to the fairly sparse dip data to be gleaned from the geological maps, and also provided valuable control points for the regional isopachs between key formation levels (Fig. 5b).
4.e. Modelling sequence of the underlying layers
We first demonstrated the feasibility of the implementation of the controls on the key Asmari surface. Then we began the implementation of a modelling sequence to construct the whole 3D model, by deducing the other stratigraphic main surfaces from the Asmari horizon. The modelling sequence is illustrated in Figure 3. Tests on different grid configuration were aimed at optimizing the mesh to minimize storage space while remaining detailed enough to honour the different control points. The final mesh was anisotropic with a larger cell size in the structure-parallel direction and a far smaller one perpendicularly; around 5 km and 300 m, respectively (best guess values). The mesh followed the trend of the structural grain, which swept around in a wide arc open to the southwest.
To deduce the lower horizons from the key Asmari horizon, regional isopach maps were constructed using two discrete datasets: firstly, the isopach values given by well data, and secondly, those provided more indirectly through the depth conversion of seismic horizon picking. To construct the depth horizons deduced from the seismic data, any variation of time picking was chosen to be translated into the isopach maps, the velocity model remaining unchanged as previously described. The construction of the isopach maps from wells required a phase of correction of the data to account for the regional dip of the surfaces between which the isopach was calculated. In particular, wells implanted on flanks of structures can give false information on thickness values between key horizons if these data are uncorrected. The same is true for the seismic data.
For each stratigraphic interval the following sequence was therefore applied. First, the apparent thickness data from wells and depth-converted seismic sections were copied to a coarse grid. These values were corrected to take into account the dip of the layers to produce a true thickness (orthogonal to the layer). This true thickness (perpendicular to the layering) was interpolated throughout the area to generate a regional isopach map. It has to be noted that owing to the large cell size, the interpolated isopach map is smoothed and does not necessarily honour the data points. A precise matching step should be later applied.
Secondly, this isopach map was used to construct the underlying surface using a layer-perpendicular shift along the direction of the surface-normal vector of the value of the true corrected isopach (concentric fold model; Galera et al. Reference Galera, Bennis, Moretti and Mallet2003). As mentioned above, the residual error for fit to well markers was then calculated at each control point. A smooth map of residual corrections was then created by interpolating the residuals regionally, which allowed a smooth fit for the well markers. It was considered that the effect of this regional application of a residual correction surface reflects the probable regional variations in seismic velocity compared with our very simple model. The remaining errors to marker fits in this case were much improved. This correction should be applied with care to ensure that the errors are not provoked by structural effects. For example, important erosional events linked to faults should be treated locally and should not be part of this process.
4.f. Modelling of diapirs
In the absence of good seismic images of the diapirs at depth, the modelling strategy for the diapirs was to combine a statistical study of the diapir diameter over the piercing diapirs present in the area and the local cropping out shape extrapolated to depth (Fig. 6 a–c). Jahani et al. (Reference Jahani, Callot, Frizon de Lamotte, Letouzey, Leturmy, Lacombe, Lavé, Verges and Roure2007) have already described the present-day morphology of salt diapirs, dividing them into six types. In the framework of the current study, only emergent diapirs not associated with the fault plane types B to E of Jahani et al. (Reference Jahani, Callot, Frizon de Lamotte, Letouzey, Leturmy, Lacombe, Lavé, Verges and Roure2007) have been considered for the study of geometry. Non-piercing diapirs of type A were not modelled, owing to the absence of good controlling data.
(a) Statistical study over the 73 piercing diapirs present in the area. Diapir surface area (in km2) versus stratigraphic column (in metres) for diapirs emerging in the indicated formation. Bars indicate the sampled values and black lozenges indicate the mean value for each formation. Mean formation thicknesses in the study area have been taken into account for the stratigraphic column. The study showed that on average, diapir diameter tended to be fairly constant in the Asmari Formation and in older series. Diameter increased dramatically above the Oligo-Miocene Asmari Formation in a characteristic ‘mushroom’ or ‘droplet’ shape with characteristic pinch and sub-salt overturn formations. This period of time corresponds to the Zagros collision and folding event. (b, c) Example of 3D modelling of diapirs. For each diapir in the model, we take into account the local shape cropping out. The modelled diapir shape is deduced from the surface extent extrapolated to depth using the above statistical study as homothetic factor.
To represent well diapir shapes, diapir diameters have been measured both parallel to and perpendicular to shortening direction for 73 diapirs. Data points were sorted according to the country rock into which the diapir is shown to be intruded at outcrop. The study showed that on average, diapir diameter tended to be smaller and fairly constant in the Asmari Formation and in older series. However, size increased dramatically above the Asmari Formation in a characteristic ‘mushroom’ shape (Fig. 6a). The parameter that tended to give the smoothest spatial distribution of diapir size turned out to be the area. Plotting diapir area against stratigraphic level provided a smoother curve than other measurements like ‘major axis’. This implies that at the scale of a diapir column, N–S-oriented horizontal compression will have been compensated, at least in part, by lateral (E–W) extrusion of the diapir to maintain surface area for a given period.
A standard diapir was therefore created using curves of the average area at each stratigraphic level. In order to represent each individual diapir surface, this curve-based diapir was then positioned at the outcrop at the correct stratigraphic level. The exact size and eccentricity of the ‘standard’ diapir was then adjusted to the values measured in outcrop for each individual diapir, and its axis rotated to coincide with the outcrop axis (Fig. 6b). Levels higher than outcrop were then removed to leave only those useful for construction of a cylindrical surface. The resulting diapir shapes were projected down to the tentative top Hormuz salt (phantom layer), and topography was applied to their upper surface (Fig. 6c).
In this model, the diapir columns were represented by a vertical shaft. This is not always supported by the analogue models (e.g. Callot et al. 2011), presenting sometimes shafts that are detached from the salt source, as well as strained, oblique or asymmetric columns, except along wrench zones. The recent analogue model results will compel us to adapt the 3D modelling to depict this behaviour in the future. The support of analogue modelling is crucial to represent the ‘column’ of the diapirs as the shape at depth is not well imaged by seismic (Callot, Jahani & Letouzey, Reference Callot, Jahani, Letouzey, Lacombe, Lavé, Verges and Roure2007; Callot et al. 2011). In offshore areas, the shaft appears to be very straight to the surface, but in our area of investigation, the seismic quality becomes generally very poor close to the salt, making precise modelling difficult.
5. Results and discussion
5.a. Uncertainty and limits of the model
The geomodel is a picture of the knowledge at a given time. It contains both ‘hard data’ such as outcrops, or well marker information, but also ‘soft data’ available in the form of depth-converted seismic lines, and even local interpretation coming from surface geology interpretations. The uncertainty comes mainly from the heterogeneity of the dataset, from the data distribution on the map (sparse wells and seismic) and from the low quality of some seismic lines. However, the geological surface dataset is of great quality, with exceptional outcrops. The geological maps are quite homogeneous, at least for the Cretaceous and Lower Tertiary formations. On the contrary, stratigraphic horizons older than Cretaceous are nearly never cropping out, and control relies only on wells and seismic interpretation. Moreover, except for two small outcrops, the Palaeozoic is not controlled by direct surface data, and its modelling relies largely on seismic interpretation.
The seismic data available only consist of sets of 2D lines, most of them quite old and showing interpretations which could be revised following any modern acquisition. In addition, except for in a small area in the extreme south of the study zone, there is no magnetic or gravimetric data of good quality, which would have helped in the modelling of the deep structures. Picking of the seismic data outside of well control allows confident determination only of a time isopach. It is impossible to determine if we need to update the models, and which of the velocity model or the isopach model should be updated, or whether we observe a joint effect (Fig. 7). In the model, it has been decided that the variations in time should be converted into thickness variations, keeping constant velocities per layer. This assumption is supported by the fact that the facies on the regional scale were quite homogeneous for the Upper Palaeozoic and the Mesozoic, being deposited on a stable carbonate platform. On the contrary, it was felt that the thickness variations could be more variable, driven by depocentres that could be provoked by the salt withdrawal (S. Jahani, unpub. Ph.D. thesis, Univ. de Cergy Pontoise, 2008; Jahani et al. Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009). In practice, the reality is probably a combination of both effects.
Illustration of the two possibilities to interpret increasing interval between seismic horizons. Left-hand side: seismic line showing increasing gap between horizons from left to right. Upper right-hand side: first modelling possibility. It assumes that real thickness of the layers is constant and that the internal velocity of the layers is decreasing to the right. Lower right-hand side: second modelling possibility. It assumes that the velocity field is constant and that the thickness of the layers is increasing to the right.
The isopach maps generated by the velocity model sometimes do not exactly honour the ‘hard data’ points since the large cell size tends to smooth the regional surfaces compared to reality. Residual errors in well marker fit ranged from 0 to 25 m, essentially owing to the difficulty of honouring a given geometry at small scale with a still relatively coarse mesh. Finer meshes will enhance the local fit.
5.b. Geological results
Building a 3D consistent geomodel requires a large investment in time. However, referencing data from various sources and from different types in three-dimensional space raises questions which would otherwise never have been apparent, and which could be critical for exploration. The modelling process allowed us to take a regional approach to the treatment of specific geological features, such as tectonic control of halokinesis, folding and basement topography.
5.b.1. Tectonic control on halokinesis
One of the main features of the studied area is the existence of numerous emerged or buried salt diapirs, composed of Upper Precambrian Hormuz salt. The studied piercing diapirs in the area display different shapes. Circular diapirs are found in the offshore area or in the coastal area, in the thrust belt front. By contrast, within the fold–thrust belt, most of the diapirs are squeezed laterally, and in the inner part only salt welds crop out.
Diapirs associated with wrench faults, tear faults or thrusts, are frequent in the central Zagros provinces (Letouzey & Sherkati, Reference Letouzey and Sherkati2004) but are less frequent in the southern Fars province (Jahani et al. Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009; Callot, Jahani & Letouzey, Reference Callot, Jahani, Letouzey, Lacombe, Lavé, Verges and Roure2007; Callot et al. 2011). One of the exceptions, the Khain Muran fold, has one of the most impressive oblique trends in the eastern Fars (Fig. 8). It is composed of three major en échelon folds, connected by a quasi continuous scar of Hormuz salt, with three major squeezed diapirs. To the south, it is thrust onto the Mishan depocentres along the fold limbs. Here in the synclines, an impressive pile of Mishan and Agha Jari deposits evidence the local increase of subsidence related to the salt escape around the former diapir, and the local increase of accommodation allowed by this escape.
Illustration of the 3D model results on the Khain anticline. (a) Three block diagrams of the model. Surfaces are from top to bottom: top Asmari, top Kazhdumi, top Hith, top Trias, top Aghar shales, top Faraghun, deep phantom horizon. (b) Vertical view (map) of the area, showing superposition of the geological surface map and the Asmari reconstructed horizon prior to erosion (in blue, cross-hatched) and position of the cross-sections (1-2-3) of the 3D view (in black). In red (darkest grey), Hormuz Salt outcrops. Note on the block diagram the modelled diapir.
The modelling emphasized the en échelon pattern of these culminations along this fold trend. It shows that the early folds, formed during the early stages of compression and located at the former salt diapir locations, initially perpendicular to the shortening direction, were later connected by several thrust faults most probably related to deeper tectonic features such as a salt ridge or basement reverse faults (see Jahani et al. (Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009) for a purely salt-related history, and Leturmy, Molinaro & Frizon de Lamotte (Reference Leturmy, Molinaro, Frizon de Lamotte, Leturmy and Robin2010) for a basement-related control).
5.b.2. Folding
The studied area in the south of the Zagros fold-and-thrust belt is one of the typical examples of a thin-skinned folded belt detached over a thick salt layer. Specific characters encountered are low tapered faulted and folded zones, and large detachment and box fold geometry due to thick salt, bending of the belt probably owing to a faster fold propagation in the centre of the salt basin where the salt is thick and the friction the lowest (Letouzey et al. Reference Letouzey, Colletta, Vially, Chermette, Jackson, Roberts and Snelson1995; allowing generalized buckling of the sedimentary cover as proposed by Mouthereau, Lacombe & Meyer, Reference Mouthereau, Lacombe and Meyer2006, Mouthereau et al. Reference Mouthereau, Tensi, Bellahsen, Lacombe, Boisgrollier and Kargar2007 and Lacombe et al. Reference Lacombe, Amrouch, Mouthereau and Dissez2007). As a general trend, we also observed an increase in the amplitude of individual folds related to the increase in shortening, similar to what was observed in the Central Zagros by Sherkati et al. (Reference Sherkati, Letouzey and Frizon de Lamotte2006). We also noticed the evolution from the poorly deformed Persian Gulf to large but relatively short along-strike detachment or box folds in the southern part, to narrower, often asymmetric anticlines in the inner part of the belt (Fig. 9, imaged on the near-top Permian depth map).
Result of the 3D model. Near-top Permian depth map over the whole study area. Depth scale increases from warm (yellow/orange) to cold (blue/purple) colours (light to dark grey). Black spots show the cropping out Hormuz salt diapirs.
Frequent out-of-sequence thrusting and numerous parasitic folds detached over intermediate décollements are observed in the flank of such structures. The Pishvar and Bavush structures, located in the central part of the belt, are good examples of such structures (Fig. 10). Both well and seismic data exist in this area, with a good calibration up to the Upper Palaeozoic. The Bavush structure is a wide asymmetric box fold. The Triassic evaporitic level is active in this structure as a minor intermediate décollement with some disharmony above and below. North of Bavush, surface anticlines are narrow and complex, and only the seismic data allow interpretation within the synclines. Figure 10 shows an unexpected deep syncline just north of Bavush. The syncline is bounded on its northern flank by a large blind thrust dipping below a narrow fault-propagation fold and the Pishvar structure.
Illustration of the 3D model results on Bavush & Pishvar anticlines. (a) Three block diagrams of the model. Surfaces are from top to bottom: top Asmari, top Kazhdumi, top Hith, top Trias, top Aghar shales, top Faraghun, deep phantom horizon. (b) Vertical view (map) of the area, showing position of the cross-sections (1-2-3) of the 3D view and position of the seismic line b002. (c) b002 line, with interpreted horizons (1 blue: Asmari; 2 green: Kazhdumi; 3 yellow: top Trias; 4 light blue: top Aghar shale; 5 orange: near-top Faraghan).
One can observe that deep horizons could be reasonably well followed in the synclines, or on the flanks of the structures, but deep horizons are not well imaged on most of the lines on the crest of the anticlines. The depth model, calibrated with wells (Bavush in that case) and regional isopachs, allows the reconstruction of the 3D view in the absence of décollement, based on the modelling of the Asmari surface.
Contrary to the central Zagros or to the north of the Qatar Arch, where thick Hormuz salt does not exist, the southern Fars area is characterized by an irregular along-strike shape of the collision-related detachment folds, with frequent apparent bending, strongly influenced by the presence of pre-existing salt diapirs (Letouzey & Sherkati Reference Letouzey and Sherkati2004; Callot, Jahani & Letouzey, Reference Callot, Jahani, Letouzey, Lacombe, Lavé, Verges and Roure2007; Jahani et al. Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009; Callot et al. 2011). A typical example of such diapir-fold influence is the Safid Herang and southern region (Fig. 11). Despite the narrowness of the syncline in between the Herang and Safid anticlines, the continuous beds at the surface do not suggest the emergence of the thrust plane between the folds. But from seismic sections the Safid anticline is a detachment fold above a thrust fault developed between two diapirs.
Illustration of the 3D model results on Safid and Herang anticlines. (a) Three block diagrams of the model. Surfaces are from top to bottom: top Asmari, top Kazhdumi, top Hith, top Trias, top Aghar shales, top Faraghun, deep phantom horizon. (b) Vertical view (map) of the area, showing position of the cross-sections (1-2-3-4) of the 3D view, and position of the seismic line U52. (c) U52 line, with interpreted horizons (1 blue: Asmari; 2 red: Kazhdumi; 3 brown: top Hith; 4 green: top Trias; 5 light blue: top Aghar shale; 6 orange: near-top Faraghan).
Similarly to the previous Bavush example, the quality of the seismic is low in the core of the anticline, but the isopach model allows building a likely 3D image of the structure.
The Gavbast anticline is one of the largest anticlines in the studied area (Fig. 12). Its rounded shape expresses the rejuvenation, during shortening, of the ascent of an underlying salt dome with a thinning of the layers. It is considered a typical example of an Early Palaeozoic buried salt dome (Jahani et al. Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009). At this spot, the huge uplift allows the observation of the stratigraphic succession down to the Fahlian beds (Lower Cretaceous). The crest of the anticline is dissected by a set of radial normal faults related to the collapse of the core of the anticline. The Gavbast dome connects three major anticlines, which are all limited along their southern flanks by reverse faults forming a zigzag structure.
Illustration of the 3D model results on the Gavbast anticline area. (a) Three block diagrams of the model. Surfaces are from top to bottom: top Asmari, top Kazhdumi, top Hith, top Trias, top Aghar shales, top Faraghun, deep phantom horizon. (b) Vertical view (map) of the area, showing superposition of the geological surface map and the Asmari reconstructed horizon prior to erosion (in blue, cross-hatched) and showing position of the cross-sections (1-2-3) of the 3D view.
5.b.3. Morphology of the pre-salt
The morphology of the pre-salt layers is one of the remaining questions which could not be solved by the available data. The top and base Hormuz salt horizons are only visible on a few offshore seismic sections, but nowhere in the onshore studied area. It is therefore not possible to identify the top Hormuz salt horizon onshore. This is probably owing to the depth of the Hormuz layer, and also to the low velocity contrast with the overlying rocks (inferred interval velocity of 4600 m/s for the Lower Palaeozoic). The deepest horizon in the model must not be interpreted as the real top Hormuz horizon, but as a reasonable minimum depth for the top Hormuz horizon. We are quite convinced that halokinesis was very active at that time and the Early Palaeozoic interval is highly anisopach (see Al-Syabi, Reference Al-Syabi2005 for the Oman salt basin and Jahani et al. Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009).
In the absence of good data concerning the thickness of the Lower Palaeozoic, the only indicators for deep basement structures is the comparison, between adjacent synclines, of the maximum depth of a given horizon. So, at a regional scale, the structures in the inner part of the folded belt are uplifted compared to the structures along the coast. This has been interpreted by several authors (see Leturmy, Molinaro & Frizon de Lamotte, Reference Leturmy, Molinaro, Frizon de Lamotte, Leturmy and Robin2010; Molinaro et al. Reference Molinaro, Leturmy, Guezou, Frizon de Lamotte and Eshragi2005 and references therein) as evidence of deep-seated basement thrusts occurring at a later stage of the Zagros orogeny in the inner part of the belt. Another regional shallowing is observed towards the west in an area where cropping out Hormuz diapirs are absent. It probably highlights the northern extension of the Qatar Arch below the thin-skin folded belt. Based on a recent gravimetric survey, a N–S-trending structural basement high is expected below the salt in the western region on the border of the Qatar Arch (Sepher, Mirhashemi & Yavari, Reference Sepher, Mirhashemi and Yavari2009). These points are in accordance with a limited but real involvement of reactivated and newly formed basement structures during folding (Mouthereau, Lacombe & Meyer, Reference Mouthereau, Lacombe and Meyer2006).
6. Synthesis and conclusions
6.a. Modelling strategy
The 3D structural modelling of the southern part of the Zagros fold-and-thrust belt represented a technical challenge. The studied area and the structural model, covering a total surface area of approximately 60000 km2, included 70 anticline structures, numerous faults and over 70 emergent Hormuz salt diapirs. Technical and geological problems faced included taking into account the numerous and closely spaced fold structures, the coexisting emergent and buried salt diapirs, along with the variable thickness and sedimentary facies of the strata, in an area with sparse well control and seismic data of generally poor quality. Fortunately, outcrop conditions are exceptional and the modelling process greatly benefited from the results of helicopter and field surveys covering most of the structures and the diapirs of the area.
For the purposes of the modelling, the data was separated into hard and soft, or interpreted, data. Hard data are not subject to review (topography, geological maps, dip and strike values measured in the field and stratigraphic marker depths obtained in the wells), whereas soft or interpreted data are subject to future changes, which will necessarily require a future update of the 3D model (2D seismic time sections and structural underground cross-sections).
For the construction of a 3D depth-domain geological model from sparse data, we had to develop a modelling strategy, both for stratigraphic surfaces and for the diapirs. For the surfaces we created a model of a key horizon for which controls could be maximized: the top Oligo-Miocene Asmari carbonates. In order to integrate the Top Asmari seismic picks, wells in the study area and neighbouring regions were used for regional velocity analysis in order to convert the seismic picks from depth to time. Horizons older than the Top Asmari were then constructed using two datasets: isopach values from well data and depth-converted time isopachs from seismic picking.
In the studied area, there are few wells and a lot of different vintages of seismic surveys. These surveys are generally not directly connected with one another and few of them are tied to wells. As a consequence, the first interpretation is often based on the seismic character. This interpretation is subject to errors owing to the possibility of mixing up the ambiguous characters of the Asmari, Khazhdumi or even Triassic seismic markers.
The 3D model is well constrained by the regional surface model deduced from the Asmari surface and well data. During the construction, the 3D interpolated surfaces could be reconverted into time, using a velocity model, and compared with previous seismic interpretations. This exercise obliged us to revise some early interpretations of seismic lines that were badly tied to wells. The 3D modelling therefore clearly improves regional interpretation. In addition, the 3D model is the only tool that allows drawing consistent cross-sections in areas where there are no seismic lines.
6.b. Geological significance
The results of the 3D geomodelling work undertaken so far in the eastern Fars domain concern seven regional stratigraphic surfaces modelled in depth, ranging from the top of the Asmari Formation (Oligocene to Lower Miocene age) to the Lower Palaeozoic above the Hormuz salt, and which respect the available hard data. Other interpretations, if carried out, must also fit these data.
The resulting geomodel merely represents the state of knowledge at a given point in the interpretation. Uncertainty arises principally from dataset heterogeneity, from the patchy seismic data distribution and poor seismic quality, and sparse well control. However, the exceptional outcrops have allowed the generation of a high-quality set of surface data, especially concerning the Cretaceous and Lower Tertiary horizons. Older horizons rely more heavily on subsurface data and uncertainty therefore increases downward in the stratigraphic column.
6.b.1. Diapir shape
Since the seismic images did not allow control on diapir shape with depth, a regional statistical study was carried out to model a standard diapir shape, which was then adapted to the individual eccentricity and stratigraphic level of each of the cropping out diapirs on a case by case basis (Fig. 6). Both field observations of the emergent salt diapirs and the statistical study performed for the purpose of the modelling, confirm some previous observations such as the progressive squeezing of the diapirs from the front to the interior of the belt, or the relations between elongation of diapirs related to thrusts or local wrench faults (Letouzey & Sherkati, Reference Letouzey and Sherkati2004; Callot, Jahani & Letouzey, Reference Callot, Jahani, Letouzey, Lacombe, Lavé, Verges and Roure2007; Jahani et al. Reference Jahani, Callot, Frizon de Lamotte, Letouzey, Leturmy, Lacombe, Lavé, Verges and Roure2007, Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009; Callot et al. 2011 and references therein).
One of the surprising results of the study was the observation of the increasing diameter of the diapirs at the time of the Zagros collision and folding event, with growth strata and overhang on the flanks of the diapirs. It shows that most of the diapirs were nearly emerged at that time, with a cone shape structure at the base of the salt (‘droplet’ or ‘mushroom’ shape of the diapirs). This event is related to a change in stress regime and sedimentation, with the arrival of clastic sediments associated with the onset of the foreland basin development as soon as the collision started. There is a competition between increased sedimentation and salt extrusion-erosion at the surface. The salt, as an incompressible medium, was squeezed upward more rapidly, changing the balance between upward movement of the diapirs and sedimentary burial. Such a correlation between an increase in salt extrusion rate and the onset of compression is also observed in other geodynamic environments associated with salt diapirs such as the West African margin or the Gulf of Mexico. Moreover, the sedimentation evolved from carbonate to clastic facies, the latter being less competent, allowing for more important salt extrusion. Nevertheless, in the absence of any data allowing the precise description of the shape of the salt diapirs at depth, it has been chosen to model the shaft of the diapir as a vertical elliptical cylinder. From analogue models in such compressive environments it is probable that the shaft is not vertical but tilted and more complex (Callot et al. 2011). This point would require further work and data.
The folded zone is characterized by large detachment folds detached above a thick Hormuz salt layer, passing to the north to out-of-sequence thrust faults. Bending of the belt is probably owing to faster fold propagation in the centre of the salt basin, where the salt is thick, and the friction the lowest (Letouzey et al. Reference Letouzey, Colletta, Vially, Chermette, Jackson, Roberts and Snelson1995). The more characteristic feature of the central part of the study is the high influence of pre-existing salt diapirs on the shape and orientation of the fold (Letouzey & Sherkati, Reference Letouzey and Sherkati2004; Jahani et al. Reference Jahani, Callot, Letouzey and Frizon de Lamotte2009; Callot, Jahani & Letouzey, Reference Callot, Jahani, Letouzey, Lacombe, Lavé, Verges and Roure2007; Callot et al. 2011).
6.b.2. Modelling of the depth to basement
Using the maps obtained from the 3D modelling (Fig. 9), we observe that there is a regional shallowing of the deeper part of the synclines immediately to the west of the studied area, probably highlighting the northern prolongation of the Qatar Arch below the thin-skin folded belt. This area presents no salt diapirism, which is another hint indicating a regional palaeo-high.
The analysis of the deep part of the synclines, based on the 3D near-top Permian map, permits us to identify some palaeo-lows and highs situated in our study area. It seems that this early structure is NE–SW oriented, with a wavelength of about 30 to 40 km.
As a matter of fact, just to the east of the Gavbendi high (prolongation of the Qatar high area), at the extreme western part of our study area, the base of the synclines are deep, indicating a palaeo-low area. This area is immediately followed by a local high (that was previously identified by Mr Sepher). Then, to the east of this zone, where diapirs are abundant (thick probable salt), there is another low area.
On the contrary, there is no observed deepening of the synclines to the north as one could expect for a fold belt zone. This absence of deepening of the near-top Permian horizon to the north may indicate an absence of flexure of the basement if we consider that there is no dramatic increase in the thickness of the Lower Palaeozoic. This is dramatically different from the expected 3° tilt of the flexed basement approaching the collision zone.
In that case, the structures of the inner part of the fault belt would be uplifted compared to the structures along the coast owing to the influence of reverse faults compensating for the flexure by stepping. The seismology indicating that there is active faulting in pre-salt layers may support this interpretation. This has been previously interpreted as deep-seated basement thrusts following the thin-skinned phase or coeval with it (Roustaei et al. Reference Roustaei, Nissen, Abassi, Gholamzadeh, Ghorashi, Tatar, Yamini- Fard, Bergman, Jackson and Parsons2010). This observation deduced from the 3D modelling is in accordance with the involvement of the basement in the shortening (Roustaei et al. Reference Roustaei, Nissen, Abassi, Gholamzadeh, Ghorashi, Tatar, Yamini- Fard, Bergman, Jackson and Parsons2010), also suggested from mechanical modelling (Mouthereau, Lacombe & Meyer, Reference Mouthereau, Lacombe and Meyer2006) or geodynamic modelling of the whole collision area (Hatzfeld & Molnar, Reference Hatzfeld and Molnar2010)
6.c. Further steps
In the scope of future work to improve the geological model, it is a definite priority to use regional understanding of the palaeogeography of the study area to refine the construction of the isopach maps, leading to a more geologically consistent regional framework. This is particularly true for the syn-kinematic post-Oligocene Asmari sediments, with regional facies changes between clastic, carbonate and evaporitic facies. Such facies changes could influence seismic velocity, and so modify locally the depth of the synclines.
The geomodel should be improved in the future by the addition of any new data such as new wells, new seismic acquisition and seismic reprocessing in order to further improve the interpretation of the structures. In such a case, the interest of a regional surface model which respects all the other available data becomes evident: converted into time, it can be used to validate and guide regional seismic interpretation, and could even be used as an a priori geological model in the case of migration in depth of the seismic processing. Additionally, the surfaces themselves can be directly used for construction of new detailed geological cross-sections of a given area, or to evaluate the size and location of potential plays for exploration, and also to model the history of burial and structuring of the region with the aim of understanding the timing of hydrocarbon generation and migration with respect to the formation of the structures.
Moreover, such thickness maps, showing the regional as well as local depocentres, will help to define the relative migration of the Hormuz salt and illustrate the relative importance of the sedimentary dynamic, the tectonic imprint and the salt withdrawal in building the sedimentary pile. This task should be undertaken jointly with a better definition of the basement geometry. In the absence of regional precise gravity and magnetic data, there is a need for new deep seismic acquisition to solve the uncertainty in the thickness of the Lower Palaeozoic and base Hormuz salt depth and morphology. But in the case of future gravimetric acquisition, the model could be used for gravimetric 3D model calibration.
Acknowledgements
The authors are grateful to the National Iranian Oil Company Exploration Directorate for supporting field trips and helicopter reconnaissance and for permission to publish this paper. We thank particularly S. M. Mohaddes, M. A. Zadehmohammadi and D. Baghbani. Thanks to NIOC geologist S. Sherkati for field work participation and discussion. We acknowledge the detailed reviews and comments of F. Mouthereau, an anonymous reviewer and editor O. Lacombe, which improved the manuscript a lot.
