3-D modeling of geologic maps from surface data
AUTHO RS
Damien Dhont $Formation de Recherche en Evolution 2639–Centre National de la Recherche Scientifique:Imagerie Ge ′ophysique,CURS-IPRA,Universite ′de Pau et des Pays de l’Adour,Avenue de l’Universite ′,Pau Cedex,France;damien.dhont@univ-pau.fr
Damien Dhont is an assistant professor of structural geology and remote sensing at the University of Pau.He received his Ph.D.in 1999from Paris 6University.His research interests include the structure of thrust and fold belts,the extensional collapse of orogens,and the study of the relationships between tectonics and volcanism.His main studied areas correspond to recent orogenic belts and plateaus of the Middle East and South America.
Pascal Luxey $Dynamic Graphics Inc.,1015Atlantic Avenue,Alameda,California,94501-1154
Pascal Luxey received his Ph.D.in 1995from Paris 6University (France).He was then granted funds to study at the Southampton Oceanog-raphy Center (United Kingdom)the interactions between the Icelandic hot spot and the Mid-Atlantic Ridge between 1995and 1998.For the last 7years,he has been a geological appli-cation specialist for Dynamic Graphics Inc.Jean Chorowicz $Unite ′Mixte de Recherche 7072–Centre National de la Recherche Scienti-fique:Laboratoire de Tectonique,case 129,Universite ′Paris 6,4place Jussieu 75252Paris Cedex 05,France
Jean Chorowicz is a professor of tectonics,remote sensing sciences,and their applications at Paris 6University.Formerly a field geolo-gist in Yugoslavia,he has concentrated on the development of remote sensing in geology since 1972.He keeps a constant research activity on orogenic belts,especially in the peri-Mediterranean region.He has studied rifts and basins,transform faults,and relationships between tectonics and volcanism.
3-D modeling of geologic maps from surface data
Damien Dhont,Pascal Luxey,and Jean Chorowicz
ABSTRACT
Recent discoveries in earth sciences are mostly related to tech-nologies allowing graphical representations of volumes.We pres-ent a way to produce mathematically and geometrically correct three-dimensional (3-D)geologic maps consisting of the volume and shape of all geologic features of a given area.The method is innovative in that it only uses surface information based on the combination of a standard geologic map,a satellite image,and a digital elevation model.It is based on a modeling algorithm that only uses surfaces calculated from scattered data points and that intersects them following a series of geologically sound rules.The major advantage of using such technology is that it provides the user with a way to quantify geology.To illustrate how a 3-D geologic map can be computed,we explain the steps taken to build a dummy model with simple faulting and depositional sequencing.The case study chosen to illustrate the method is the Beirut watershed (Lebanon),an area with relatively simple geology.The 3-D visualization and cross sections help in the understanding of the geometrical relationship between the different geologic fea-tures,allowing a reexamination of the tectonic history of the area during the late Mesozoic.
INTRODUCTION
Modern geology requires accurate representation of layer volumes.Three-dimensional (3-D)geologic models are increasingly the best method to constrain geology at depth.Until now,geologic volume modeling has been based on the interpretation of expensive two-dimensional (2-D)and 3-D seismic surveys and/or well-log data.Hence,it has been typically used by the oil industry for explo-ration and production.
GEOHORIZONS
AAPG Bulletin,v.89,no.11(November 2005),pp.1465–14741465
Copyright #2005.The American Association of Petroleum Geologists.All rights reserved.
Manuscript received October 12,2004;provisional acceptance January 12,2005;revised manuscript received June 6,2005;final acceptance June 27,2005.DOI:10.1306/06270504108
Quantitative structural analysis related to surface geologic map-ping has been obtained using stereoscopic remote sensing imagery (Berger et al.,1992;Bilotti et al.,2000).It allows the calculation of the strike and dip of bedding,which can be used to supple-ment field mapping (Banerjee and Mitra,2004).In a recent ar-ticle,Fernandez et al.(2004)presented a methodology to generate 3-D reconstructions of geologic surfaces by integrating geologic mapping and field data.Here,we go further in that we present a way to generate efficient and accurate representations of geologic structure volumes at depth using existing and inexpensive data sets.This approach is innovative in that it is based only on surface information from published geologic maps,remote sensing data,and a digital elevation model (DEM).
The method is applied to the Beirut watershed (Lebanon),which is characterized by well-exposed and well-mapped structures (Figure 1).The area is located on the western side of the Lebanon ridge and is underlain by Late Jurassic to early Cenozoic interbedded sandstone and fissure-karstic limestone that are cut by volcanic intrusions.All layers are offset by subvertical faults,with small amounts of throw,striking mainly northeast and east-southeast (Freund et al.,1970).East of the watershed,the strata are almost flat.Near the coast,they abruptly dip west.This peculiar shape,called the western Lebanon flexure,is dated to at least early Mio-cene (Walley,1998).
Our aim is to present a 3-D geologic map that gives quan-titative data not only at the surface (i.e.,bed strikes and dips)but also at depth.The map then provides the opportunity to describe the geometry of the geologic structures through different graphical representations.
STATE OF THE ART
The concept of a volumetric or 3-D geologic map needs to be clarified because it is commonly misused.Fernandez et al.(2004)present a good review of the different methods developed for the generation of 3-D reconstructions.On a standard paper geologic map,geologic features such as lithologic boundaries or faults rep-resent four variables:the spatial coordinates (x ,y ,z ),expressed on the topographic map,and another variable that represents what-ever criteria are used by the geologist to differentiate geologic units (time,facies,depositional event,etc).Hence,the geologic map corresponds to a planar view (2-D)representing four-dimensional variables.Digitizing geologic boundaries on a map creates a 2-D matrix (x ,y )containing elements corresponding to an isochron (t ),but the elevation (z )is not considered.Draping the scanned geo-logic map onto a DEM is a way to account for the missing variable (z ).Such technology allows a better representation of the geology.Because a 2-D topographic map represents the variation of (z )along the x -and y -axes,a flying-carpetlike view of a DEM remains a 2-D object despite the fact that it is seen in 3-D.Displays of
ACKNOWLEDGEMENTS
This study has been made in the frame of the collaboration between the Saint-Joseph Uni-versity in Lebanon and the University of Pau in France.This project has been funded by the Research Council of Saint-Joseph University.The article benefited from constructive and thorough reviews by Jim Granath,Sandro Serra,and an anonymous referee.The English has been improved thanks to careful corrections of Carol Man.
1466
Geohorizons
scanned geologic maps draped onto DEMs in rotating perspective views have therefore been misleadingly called 3-D geologic maps.Such representations have allowed more information to be gathered on maps,but this still cannot provide any volumetric underground information as a true 3-D geologic map should.
The advent of geographic information systems (GIS)represented a major step forward in the reali-zation of maps by providing the ability to include and superimpose several data types.This tool,how-ever,only allows the manipulation of 2-D maps,which,when superimposed,do not produce a clear image of the subsurface,as would be the case for a 3-D map.
A 3-D map is the representation of geologic units and structures seen as actual volumes following a
3-D
Figure 1.Simplified geologic map of Lebanon and location of the studied area (dark frame)(modified from Dubertret,1955,1975,in Walley,1998).
Dhont et al.
1467
matrix(x,y,z)at a given time(t).A geologic unit vol-ume is defined within two surfaces corresponding to the layer’s top and bottom.The volume can be further limited either by fault surfaces(if the geologic unit has been faulted)or by erosion surfaces.Although such maps can be constructed from2-D and/or3-D seismic surveys,they are difficult to acquire in areas of rough topography.The originality of this approach to digital mapping is that a true3-D map is built only from surface data extrapolated to depth.A true3-D map represents the only way to precisely quantify geology by providing,for example,volumetrics and enveloping surfaces.
Our new technique takes full advantage of a flexible package(EarthVision1)combined with a data set that includes the definition of strike and dip of horizons and faults.This work is unique in that it al-lows the accurate modeling and volumetric represen-tation of any structural features in any geologic setting. The slightest geologic detail,namely,fault offsets or changes in dip directions,captured from the input data set can be rendered with accuracy in the3-D map.The faithfulness of the final display demonstrates that no simplification needs to be applied to the input data.The fact that all intersections between horizons and faults are automatically calculated renders the study feasible in a bearable time frame and with a possibility to up-date the model very easily.Before using this package and this technique,intersections between faults and horizons had to be designed grid node by grid node, making the study of complex faulting impossible;the fault hierarchy was simply too complex to define in3-D. Other packages simply cannot render reverse faults or Y-shape faults.In these cases,fault verticalization is the only option.The rendering of such models as volumet-ric geologic maps,however,leads to errors and mis-interpretations.Areas with very simple geologic set-tings could only be modeled,or drastic simplifications of complex areas had to be applied,which defeats the purpose of such studies.Moreover,the surface ren-dering of other packages does not allow the clear un-derstanding that our volumetric map provides. MODELING TECHNIQUE
To illustrate the modeling technique,we built a simplistic dummy model based on a geologic map (Figure2a).The geology of the model corresponds to the superimposition of three layers(yellow,red,and green)all cut by a single,nearly vertical fault.To gen-erate a DEM,we first digitized the topographic con-tours from the geologic map.They are represented by the black and green lines(Figure2a)and small green dots(Figure2b).We then calculated the DEM from these points.The elevation value of each digi-tized geologic feature contour is calculated by its pro-jection onto the DEM in what is called a‘‘back-interpolation’’process,i.e.,the interpolation of the DEM elevation at the(x,y)locations of the points forming the digitized contours.The fault surface (i.e.,2-D grid)is defined first and used to split the 3-D model space into fault blocks(Figure2c).Using a cluster of points in each fault block,a numerically defined surface is then passed through these points to predict the position of each geologic feature through-out the fault block volume(Figure2d).The red ho-rizon is represented in solid color in one fault block and transparent in the other.Its volume is limited by three different surface groups:the model limits,the horizon limits(top,bottom,and DEM),and the fault surface.We used the3-D modeler to intersect the 2-D grids following geologic rules and to generate layer volumes.The whole volume is assembled from the bottom up.The oldest layer is modeled first,its volume defined by its top,eventually by the fault surfaces when present and by the base of the model. The overlying layer is defined in the same way,and
Figure2.Steps used to extrapolate a2-D geologic map into a3-D map.The input data(a)are the geologic map alone.It depicts three layers cut by a single,nearly vertical fault.The black crosses show the points digitized for each geologic feature:the layer tops and the fault trace.The black and green elevation contours have been digitized to generate the DEM.All data used on the model are shown in(b).A back-interpolation onto the DEM is applied to define an elevation to each digitized contour represented by yellow, red,and green dots.Digitized fault points are represented in pink,and the corresponding2-D grid is shown in(c).The fault splits the model into two fault blocks.The red horizon is represented in solid color in one fault block and transparent in the other(d).Its volume is limited by different surfaces:the model limits,the top of the red layer,the top of the yellow layer,the fault surface,and the DEM.The horizon volumes are defined by their top surfaces;the top of the underlying horizon is therefore the base of the horizon above(e).The resulting model(f)shows the relationships between the layers and the DEM because this latter surface has been used as an erosion to truncate all underlying volumes.
1468Geohorizons
the base of its volume consists of the top of the un-derlying one (Figure 2e).This process is applied for the entire geologic sequence.The DEM truncates each underlying volume and fills all remaining spaces below it,giving the model its 3-D volumetric map shape (Figure 2f).
Once the preliminary settings are defined in the package (mainly geologic rules,fault surfaces
shape,
Dhont et al.1469
and fault hierarchy),the computing of the model takes only minutes.The validation of the model is done in a trial-and-error process until each geologic layer top fits its contour on the geologic map,and the global geometry makes geologic sense.Modifications consist of moving,adding,or removing input scattered data points.The back-interpolation of the digitized geo-logic contours onto the DEM is noisy(abrupt changes in dip and dip direction)when the contours are defined by too many points.We experienced that a few well-picked points giving both a clear dip and location(for example,the valley-shaped contour in Figure2a)lead to a better defined surface.
CASE STUDY
Workflow
The available data used to model the Beirut watershed are inexpensive or even freely downloaded via the In-ternet(Figure3).They consist of(1)a30-m(100-ft) ground pixel DEM,derived from the Terra spacecraft Aster instruments(EOS Web site,http://edcimswww https://www.wendangku.net/doc/a115859637.html,/pub/imswelcome/);(2)a color composite Landsat5Thematic Mapper image at28.5-m(93.5-ft) pixel size(NASA Web site,https://www.wendangku.net/doc/a115859637.html, /mrsid);and(3)the scanned1:50,000Beirut geologic map(Dubertret,1951).The georeferencing of both the geologic maps has been done using ground control points referenced from the already georeferenced DEM.We picked three ground control points evenly spread on the study area leading to an error less than 100m(330ft)in the(x,y)direction.The error in elevation is directly linked to the quality of the DEM and lies within a few tens of meters.Considering the size of the studied area,these values are acceptable.A more detailed study would have required better input data,such as a more accurate DEM and a higher reso-lution image.The power and uniqueness of this real case study lie in the fact that the3-D model is fast to build and easy to update.The first-pass model took less than20hr to set up and compute.Its fine tuning took about the same time.Updating it takes only a few minutes because the whole model structure is already defined.In this process,the geologic surfaces are de-fined solely from their intersection with the topo-graphic surface.The layer geometries at depth, therefore,are obviously extrapolated and may not correspond to the actual geology.Hence,a model built this way is only valid close to the surface.Further im-provements can be made if field measurements,well surveys,cross section interpretations,etc.are brought into the modeling process.These measurements pro-vide valuable information that can be easily added to the model.
The finalized model was validated once there was a close to perfect match between the geologic map of Dubertret(1951),the satellite image,and the inter-sections of the calculated geologic surfaces with the DEM.The technique used to build the Beirut water-shed model follows exactly the steps described in the section on Modeling Technique.The main difference resides in the combined use of the satellite image and the digitized geologic map.The contours were first digi-tized from the geologic map and then corrected follow-ing the geomorphic interpretation from the image draped onto the DEM.The rock contrasts are well highlighted on the image,thus allowing a more precise drawing of the geologic contours than from the map alone.One difficulty resided in the modeling of fault surfaces.The fault traces on the map are mostly rep-resented by straight lines leading to the modeling of vertical surfaces.We corrected some of the fault bed-dings according to the geologic setting by adding points to shape their surface.Once the fault geometry and hierarchy is set,the layer contours have to be adjusted (1)at depth,because away from the digitized con-tours,the layer tops are extrapolated by the package;
Figure3.Data used for the construction of the3-D geologic model.The topographic surface is derived from the30-m(100-ft)DEM obtained from Terra Aster Product,freely available on the Internet.3?vertical exaggeration.The digitized geologic contours of the 1:50,000scale sheet of Beirut(Dubertret,1951)are corrected and georeferenced using the DEM in an interactive process.The corresponding2-D grids are calculated and intersected following geologic rules using an earth-modeling package.Each different geologic unit is then represented as a volume.A layer is made of an assemblage of several fault block volumes.The perspective view represented here corresponds to the digital geologic contours and faults.Two cross sections are shown and display the relationship between horizons and faults.Faults are nearly vertical,mostly with a normal component as is depicted on the top of the horizon J5 corresponding to the basement of our geologic sequence for this model.The depth color contouring highlights the deformation of the J5horizon along faults.
1470Geohorizons
Dhont et al.1471
1472Geohorizons
and(2)at surface level,because inaccurate point lo-cation leads to erroneous local bedding.
In the western part of the model,the distribution of points was sufficient to prevent the need for any adjustments;therefore,the thickness changes and cur-vatures of the surfaces are strictly caused by the local dip angle and dip direction calculated from the contour themselves.Some layers with limited outcrops were modeled as intermediate surfaces taking their shape from another layer defined by a better set of points. This is the case for the J7(Tithonian)layer that out-crops only in the northeast of the model,whose shape is derived from the J6(Oxfordian–Kimmeridgian). In the eastern part of the model,the shapes of the surface below the sea have been modeled in the sim-plest way possible.The extrapolated part of the ho-rizon surfaces has been adjusted by adding a few points at depth below the sea to constrain the dip of all sur-faces to remain constant and to give the horizons a constant thickness.
Analysis
The3-D map makes the geometry of the area under-standable at first glance,even by nonspecialists,the relationships between layers and faults being self-explanatory(Figure4).Such a3-D map displays in-formation on the regional tectonic history,and it is provided in a way that a simple standard geologic map cannot display.Graphical representations allow the user to examine it from various directions,slice it to generate cross sections,or disassemble it to examine individual geologic units.We have cut the3-D model into four pieces along two cross sections.On the east-west cross section,the C1(Valanginian–Hauterivian) layer progressively thickens from west to east because of the subsidence of the eastern part of the area.The 3-D geologic model shows numerous,nearly vertical faults with a couple of hundred meters throw,bound-ing grabens a few kilometers in size.The thicknesses of the C2a(Barremian)and C2b(Aptian)units vary from one fault block to another,recording the synsed-imentary fault motion and the filling of the small grabens during this middle Cretaceous extensional event.Syntectonic deposition during this time interval has also been described elsewhere in Lebanon and Syria (Guiraud and Bosworth,1997;Fleury et al.,1999; Sawaf et al.,2001;Homberg et al.,2003).For illus-trative purposes,we have grouped all layers younger than C2b together floating above the syntectonic units. Note that the C3(Albian)and C4(Cenomanian)ho-rizons,although discontinuous on the model caused by erosion,show a much more constant thickness,sug-gesting the end of this middle Cretaceous extensional event.
CONCLUSIONS
From a DEM,a satellite image,and a geologic map alone,we are able to generate an accurate3-D geologic map of the Beirut watershed area(Lebanon)that high-lights with extreme clarity the structure of the studied https://www.wendangku.net/doc/a115859637.html,ing these data,it is possible to quantify,with as much accuracy as the input data permit,the geomet-rical relationships between faults and horizons.The3-D geologic map reveals features that standard2-D geo-logic maps or satellite images by themselves could not show.For instance,in the studied area,we display ex-tensional deformations during the middle Cretaceous that are expressed by fault throw variations,layer thick-ness changes,and stratigraphic relationships.
This methodology can be applied to other parts of the world and to other geologic tasks,using the same kind of input geologic map,complemented with field and remotely sensed data where needed and available. Three-dimensional maps will also prove to be power-ful tools in other disciplines such as agricultural ir-rigation,well management,natural risks management, and mining.The individual units can be augmented with other attributes,such as rock properties,making it advantageous to integrate the3-D map with,for example,petrophysical modeling.
Figure4.The3-D map has been cut along two cross sections(displayed in red in the inset).Each side of the cross section is pushed away from the other to display the cross section in the lowest layers of the model and its relationship with the neighboring areas.The geologic sequence has been disassembled into two units.The top unit is shown floating and contains the layers younger than the C2b horizon.The unit below contains the remaining sequence corresponding to layers older than C3.The advantage of such a3-D map is that the geometric relationships between the geologic features are self-explanatory and easy to interpret.The west-east cross section shows(1)variations in the sediment thickness along the entire cross section for C1,(2)synsedimentary variations in the thickness of C2a and C2b from one fault block to another,followed by(3)the uniform thickness of C3and C4.
Dhont et al.1473
REFERENCES CITED
Banerjee,S.,and S.Mitra,2004,Remote surface mapping using orthophotos and geologic maps draped over digital elevation models:Application to the Sheep Mountain anticline,Wyo-ming:AAPG Bulletin,v.88,p.1227–1237.
Berger,Z.,T.H.Lee Willams,and D.W.Anderson,1992,Geologic stereo mapping of geologic structures with SPOT satellite data: AAPG Bulletin,v.76,p.101–120.
Bilotti,F.,J.H.Shaw,and P.A.Brennan,2000,Quantitative struc-tural analysis with stereoscopic remote sensing imagery:AAPG Bulletin,v.84,p.727–740.
Dubertret,L.,1951,Carte ge′ologique de Liban au1:50,000,Bey-routh sheet,Paris,Institut Ge′ographique National. Dubertret,L.,1955,Carte ge′ologique du Liban au1:200,000, Beyrouth,Paris,Institut Ge′ographique National. Dubertret,L.,1975,Introduction a`la carte ge′ologique au1/50,000 du Liban:Notes et Me′moires du Moyen–Orient,v.13,p.345–403.
Fernandez,O.,J.A.Mun?oz,P.Arbue′s,O.Favilene,and M.Marzo, 2004,Three-dimensional reconstruction of geological surfaces:An example of growth strata and turbidite systems from the Ainsa basin(Pyrenees,Spain):AAPG Bulletin,v.88,p.1049–1068.Fleury J.,J.Chorowicz,and J.Somma,1999,Discussion on trans-current fault activity on the Dead Sea transform in Lebanon and its implications for plate tectonics and seismic hazards: Journal of the Geologic Society(London),v.154,p.757–760.
Freund,R.,Z.Garfunkel,I.Zak,M.Goldberg,T.Weissbrod,and
B.Derin,1970,The shear along the Dead Sea rift:Philosoph-
ical Transactions of the Royal Society of London,v.267, p.107–130.
Guiraud,R.,and W.Bosworth,1997,Senonian basin inversion and rejuvenation of rifting in Africa and Arabia:Synthesis and implications to plate-scale tectonics:Tectonophysics,v.282, p.39–82.
Homberg,C.,E.Barrier,M.Mroueh,W.Hamdan,F.Hijazi,and C.
Mu¨ller,2003,Early Cretaceous and Eocene extensions in Lebanon:European Geophysical Society Meeting Abstracts, Nice,v.5,p.10,461.
Sawaf,T.,G. E.Brew,R.K.Litak,and M.Barazangi,2001, Geologic evolution of the intraplate Palmyride basin and Euphrates fault system,Syria:Peri-Tethys Memoir6,p.441–467.
Walley,C.D.,1998,Some outstanding issues in the geology of Lebanon and their importance in the tectonic evolution of the Levantine region:Tectonophysics,v.298,p.37–62.
1474Geohorizons