3D object-oriented image analysis in 3D geophysical modelling: Analysing the central part of the East African Rift System

(1)

ContentslistsavailableatScienceDirect

International

Journal

of

Applied

Earth

Observation

and

Geoinformation

j ou rn a l h o m epa g e :w w w . e l s e v i e r . c o m / l o c a t e /j a g

3D

object-oriented

image

analysis

in

3D

geophysical

modelling:

Analysing

the

central

part

of

the

East

African

Rift

System

I. Fadel

a,b,∗

_,

_M.

_van

_der

_Meijde

a

_,

_N.

_Kerle

a

_,

_N.

_Lauritsen

c

a_University_of_Twente,_Faculty_of_{Geo-information}_Science_and_Earth_Observation_(ITC),_P.O._Box_6,₇₅₀₀_AA_Enschede,_The_Netherlands b_Geology_Department,_Faculty_of_Science,_Helwan_University,_Ain_Helwan,_Egypt

c_DTU_Space,_National_Space_Institute,_Elektrovej,_Building_327+328_and_Ørsteds_Plads,_Building_348,_DK-2800_Kgs._Lyngby,_Denmark

a

r

t

i

c

l

e

i

n

f

o

Articlehistory:

Received15August2013 Accepted8November2013 Availableonline8December2013 Keywords:

Satellitegravity 3Dgravitymodel

Object-orientedimageanalysis Inversion

TanzaniaCraton Seismictomography

a

b

s

t

r

a

c

t

Non-uniquenessofsatellitegravityinterpretationhastraditionallybeenreducedbyusingapriori infor-mationfromseismictomographymodels.Thisreductioninthenon-uniquenesshasbeenbasedon velocity–densityconversionformulasoruserinterpretationofthe3Dsubsurfacestructures(objects) basedontheseismictomographymodelsandthenforwardmodellingtheseobjects.However,thisform ofobject-basedapproachhasbeendonewithoutastandardizedmethodologyonhowtoextractthe sub-surfacestructuresfromthe3Dmodels.Inthisresearch,a3Dobject-orientedimageanalysis(3DOOA) approachwasimplementedtoextractthe3Dsubsurfacestructuresfromgeophysicaldata.Theapproach wasappliedona3DshearwaveseismictomographymodelofthecentralpartoftheEastAfricanRift System.Subsequently,theextracted3Dobjectsfromthetomographymodelwerereconstructedinthe 3DinteractivemodellingenvironmentIGMAS+,andtheirdensitycontrastvalueswerecalculatedusing anobject-basedinversiontechniquetocalculatetheforwardsignaloftheobjectsandcompareitwiththe measuredsatellitegravity.Thus,anewobject-basedapproachwasimplementedtointerpretandextract the3Dsubsurfaceobjectsfrom3Dgeophysicaldata.Wealsointroduceanewapproachtoconstrainthe interpretationofthesatellitegravitymeasurementsthatcanbeappliedusingany3Dgeophysicalmodel. ©2013ElsevierB.V.Allrightsreserved.

1. Introduction

Invertinggravitydatahastraditionallybeenoneofthemain

geophysicalchallenges.Becauseofitsmonopolepotentialﬁeld,and

resultingmeasuredbulkproperty,retrievingaccurateobject

infor-mationis oftena matteroftrade-offsand bestguesses.Usinga

priorigeologicaland/orseismologicaldataisawayofreducingthis

non-uniquenessofgravityinversions.

Therearetwomaingenericapproacheshavebeenfollowedin

thepast.Theﬁrstistheconversionoftheseismicvelocitiesinto

densitiesusingempiricalstandardvelocity–densityrelationships

(Birch,1961;Gardneretal.,1974;ChristensenandMooney,1995; Godfreyetal.,1997;Brocher,2005)orincludingtemperature,

pres-sure,andmineralcompositioneffects(Ravatetal.,1999;Korenaga

etal.,2001;Cammaranoetal.,2011).Calculateddensitieswerethen

usedtoproducetheestimatedgravitysignalandcompareitwith

themeasuredgravity,ortoconstrainanyfurtherinversion

pro-cess(BezadaandZelt,2011).Thisapproachusuallyhastodealwith

∗ Correspondingauthor.

E-mailaddresses:i.e.a.m.fadel@utwente.nl,islam.geophysics@gmail.com

(I.Fadel),m.vandermeijde@utwente.nl(M.vanderMeijde),n.kerle@utwente.nl

(N.Kerle),nlbla@space.dtu.dk(N.Lauritsen).

uncertaintyintheseismicdataandthevelocity–densityconversion

relationshipsthatresultfromthedependencyonanduncertainty

inphysicalparameterssuchastemperature,pressureand

petro-logy(BezadaandZelt,2011;Cammaranoetal.,2011).Thesecond

approachisbasedontheinterpretationoftheseismicmodelsto

deriveobjects(GötzeandLahmeyer,1988)suchaslayers,plate

boundaries,subductionzones,anddomes.Subsequentlythe

densi-tiesoftheseobjectsarecalculatedusingoneofthevelocity–density

relationshipsorthroughatrialanderrorapproach.Afterthat,the

forwardsignalassociatedwiththeseobjectsiscomparedwiththe

measuredgravity(e.g.,Sanchezetal.,2011).Afollowingstepcan

beapplyinginversionapproachestoreducethemisﬁtbetweenthe

calculatedandthemeasuredgravity(e.g.,Ebbingetal.,2001).This

gravitymodellingobject-basedapproachhasbeendonewithouta

standardizedmethodologyonhowtoextractthesubsurface

struc-turesfrom3Dsubsurfacedata,andhastraditionallymainlybeen

basedontheuserexperienceorpre-existingdata(e.g.,Ebbingetal.,

2001;Sanchezetal.,2011).

Inthisresearchwereportonacasestudyfromthecentralpart

oftheEastAfricanRiftSystemandhowanadvanced3D

object-orientedImage Analysis (3DOOA) approach can provide more

objectivedeﬁnitionsofobjectboundariesforgravitymodellingand

inversion.3DOOAwillbeusedtoextract3Dobjectsfroma

tomo-graphyvoxelmodel.Then,theywillbereconstructedinthe3D

(2)

I.Fadeletal./InternationalJournalofAppliedEarthObservationandGeoinformation35(2015)44–53 45

object-basedinteractivegravitymodellingenvironmentIGMAS+.

Subsequentlytheirdensitycontrastvalueswillbecalculatedusing

anobject-basedinversion.Finally,thegravityresponseoftheﬁnal

modelwillbecomparedwiththesatellitegravitysignal.

2. Dataset

2.1. 3Dtomographymodel

ThestudyareaisthecentralpartoftheEastAfricanRiftSystem.

Adamsetal.(2012)developeda3Dshearwaveseismic

tomogra-phymodelofthearea,whichisvoxel-basedandcanbedividedinto

twomainsections.Theﬁrstisthecrustalpartofthemodel,which

extendsfrom0to40kmwith1kmincrements.Itrepresentsthe

crustalstructureofthemodeldowntoaﬂatMohoat40kmdepth.

Thesecondistheuppermantlepartofthemodelrangesfroma

depthof40kmdownto400kmat10kmincrements.Thispartof

themodelcanbedividedintotwomainzones,shallowanddeep,

basedontheirrelativepositioninthemodel.Theshallowzoneof

theuppermantlepart(from40kmdownto∼200km)is

charac-terizedbyhighvelocitiesoftheTanzanianCratonandUgandan

basementcomplex,andlowvelocitiesofthesurroundingtectonic

settingssuchasriftbranches, theProterozoic mobilebeltsand

Cenozoicvolcanism.Thedeepzoneoftheuppermantlepart(from

∼200kmdownto400km)ischaracterizedbydominantlow

veloc-itiesrelatedtoaheadplumeorpartofthedeepmantlesuperplume

inwesternAfrica(Adamsetal.,2012).

2.2. Satellitegravitydata

ThegravitydatausedinthisstudyareGOCEonlymodels(see

vanderMeijdeetal.(2013)fordetails)andthehighresolution

EIGEN-6C(Förste et al., 2011), a combinedgravity ﬁeld model

of LAGEOS/GRACE,GOCE and DTU10 upto a maximumdegree

andorderof1420.Thefreeairgravitydatawereobtainedfrom

the InternationalCentre for Global Earth Models (ICGEM) web

site(http://icgem.gfz-potsdam.de/ICGEM/)witharesolutionof0.1

degree(Fig.1A).Theeffectsofthenearbyterrainwereremoved

usingETOPO1topographydata(Fig.1B),alsodownloadedthrough

theICGEMwebsitewithsimilarresolution(0.1degree)tothe

mod-els.The effect of the terrain correctionon the satellitegravity

measurementswassmallwithamaximumvalueof3.16×10−5

m/s2₍_Fig.₁_C),_in_line_with_ﬁndings_by_Mishra_et_al.₍₂₀₁₂₎_._Then,

thedatawereBouguercorrected,removingtheeffectduetothe

massexcessordeﬁciencyasaresultofhighorlowelevation

rela-tivetothemeansealevel.TheBouguercorrectionalsousedETOPO1

topographydatawiththesamespatialresolution(0.1degree)asthe

satellitegravitydataandareferencedensityof2.67t/m3₍_Fig.₁_D).

3. 3DOOA

Objectorientedimageanalysis(OOA),alsodeﬁned asobject

basedimageanalysis(OBIA)orgeographicobjectbasedimage

anal-ysis(GEOBIA),isanimageanalysistechniquethatisbasedonthe

analysisofsegments,orobjects,insteadofpixels.Ingeneral,OOA

isa2-stepprocedurestartingwithimagesegmentationthatis

fol-lowedbytheclassiﬁcationofthesesegments.Thesegmentation

stepcomprisesofdividingimagepixelsintocontiguoussegments

orobjects,basedonpre-deﬁnedknowledgebasedorhomogeneity

criteria,suchascolourorintensity.Thesubsequentclassiﬁcation

stepusesthesegmentsassociatedattributeswhichcanbespectral,

geometric,textural,contextual,andalsogeologicalandphysical,to

differentiatebetweenthedifferentobjectspresentintheimages.

Consequently,OOAmimicscognitivevisualimageanalysis

tech-niquesthroughaformofknowledge-drivenanalysis(Marthaetal.,

2011).Recently,3DOOAwasdevelopedinthebiomedicalﬁeldto

extract3Dobjectsfromimagestacks(Schoenmeyeretal.,2006).

Theﬁrstapplicationto3Dsyntheticandseismologicalmodelswas

recently carried out (Fadel et al., 2013).Thatstudy speciﬁcally

appliedtheapproachineCognitionsoftwareforderivingobjects

ina3Dseismictomographymodel.

eCognitionwastheﬁrstOOAtoolwhenintroducedin2000and

sincethen about50–55%of allOOApublishedscientiﬁc studies

haveusedit(Blaschke,2010).Thenthenamewaschangedinto

Deﬁniensin2007.Afterthatthepackagewassplitintotwo

soft-warepackages;Deﬁnienssoftwarethat isrelatedtobiomedical

applications,andeCognitionthatismainlyappliedforgeospatial

applications.eCognition isbuiltonthe conceptof rulesetsthat

deﬁnetheobjectextractionprocessinwhichacomplexprocessing

schemeof segmentation/classiﬁcationalgorithmscanbe

imple-mented(Fadeletal.,2013).Therulesetisasemiautomaticsince

itisbuiltonsegmentationandclassiﬁcationalgorithms,whichare

originallydevelopedintheeCognitionpackage,thatareusedto

seg-mentortoclassifythe3Dobjectsfromtheimagestack.Theuser

onlyneedstodeﬁnethethresholdsthataresuitabletodeﬁnethe

targetsegmentationschemeorthetargetobjecttobeextracted.

3DOOAwasappliedonthe3Dtomographicmodelin2stages.

First,therulesetofFadeletal.(2013)wasusedtoextracttheobjects

fortheuppermantlepartofthemodelbetween>40–400km(Fig.2

UpperMantlePart).Theuppermantleobjectscanbeseparatedin

twozonesbasedontheirrelativedepthintheuppermantlepart

ofthemodel:

(1)The shallow zone (>40-∼200km) consists of high velocity

objects (craton (C)and shallowhighvelocity object (SHV)),

lowvelocityobject(rifts(R)),andboundaryzones

surround-ingbothhighandlowvelocityobjects(boundaryshallowhigh

velocity(BSHV)andaboundaryshallowlowvelocity(BSLV),

respectively).

(2)Thedeepzone(∼200–400km)isdominatedbyalowvelocity

object(LV)withitsinnerpart(ILV)followedbyahigh

veloc-ityobjectatthedeepestpartofthemodel(deephighvelocity

(DHV))withaboundarybetweenitandtheupperlowvelocity

object(boundarydeeplowvelocity(BDLV)).

Second,anewrulesetwasdevelopedforthecrustalpartofthe

modelfromthesurfacedowntoadepthof40km.Thenewcrustal

partofthemodel(AubreyaAdams,personalcommunication,April

29,2013)wasvoxelbasedwith1kmdepthinterval.

Thenewrulesetforthispartwasbasedonvisualinterpretation

andanalysisofthe3Dhistogram,summarizingthe3Ddistribution

ofdensityvalueswithdepth(Fig.3B),similartoFadeletal.(2013).

Thetomographicmodelhadahardboundaryof2discrete1km

layersaboveandbelowtheMoho(Adamsetal.,2012),withamajor

velocityjumpoverthesetwolayers.Therulesetmadeuseofthis

hardboundarytodistinguishthecrustalpartfromtheuppermantle

part.Thenewrulesetconsistedof2mainsteps(seesupplementary

dataformoredetails).

Intheﬁrststepmulti-thresholdsegmentationwasusedto

sep-aratethecrustalfromtheuppermantle partofthemodel.This

algorithmsplitstheimageintosegmentsandclassiﬁesthembased

onpre-deﬁnedthresholds.40kmdepthwasusedasathresholdto

separatethecrustalfromtheuppermantlepart(Fig.3A).The

subse-quentfocuswasonthecrustalpartofthemodel(0–40km).Aseries

ofassignclassalgorithmswasappliedtosplititintoﬁveobjects

usingthresholdsthatweredeﬁnedbasedonthevisual

interpreta-tionsupportedbytheanalysisofthe3Dhistogram(Fig.3B).These

5units(Fig.2CrustalPart)donotnecessarilyhaveageo/physical

meaningintermsoftheirvelocities,due tothesharpboundary

withlargevelocitycontrastattheﬁxedlayerat40kmdepthinthe

tomographymodel.However,theyhighlighttherelativevelocity

(3)

Fig.1. Thecorrectionofthesatellitegravitysignal.(A)ThefreeairEIGEN-6cdata.(B)ETOPO1topographymap.(C)TheterraincorrectionmapusingETOPO1topography andreferencedensity2.67t/m3_._(D)_The_Bouguer_anomaly_map_after_the_Bouguer_and_terrain_corrections.

valuesupto3.2km/s,andwassurroundedbyunit2which

repre-sentstheboundaryforunit1andhasavelocityrangefrom3.2km/s

up to3.4km/s. Unit 3, with a velocity range of 3.4–3.55km/s,

wascharacterizedbylowvariationinthedensitydistributionas

canbeseeninthe3Dhistogram.Unit4,withavelocityrangeof

3.55–4.1km/s,surroundedtheﬁnalunit5thatwascharacterized

byvelocitieshigherthan4.1km/stherebyrepresentingthehighest

velocityinthecrustalpart.3DOOAbehavedasanysimple

isosur-faceduetothelackofaprioriinformationaboutthecrustalpart;in

additiontothelimitedabilityofseismictomographytechniquesin

suchshallowdepthswhereotherapproachesaremorepreferable

suchasreceiverfunctions(LiuandGu,2012).

4. Forwardmodelling

4.1. Objectsreconstruction

The objects extracted in the OOA part were re-constructed

in IGMAS+ (3D Interactive Gravity and Magnetic Application

System;GötzeandLahmeyer,1988;Schmidtetal.,2010),which

facilitatesinteractiveforwardmodellingof thesubsurface inan

object-orientedenvironment.TheresultsfromOOAwereexported

asclassiﬁedcrosssectionsthatwereusedtoconstraintheshape

oftheobjectsthroughthereconstructionprocess.27sectionswere

createdinIGMAS+with50kmaveragedistancebetweenthemto

ﬁttheoriginal3Dtomographymodelthathadahorizontal

resolu-tionof0.5degree.Then,OOAclassiﬁedsectionswereimportedin

IGMAS+andco-locatedwiththecreatedIGMAS+sectionstostart

re-constructingtheobjects.Theuppermantlepartofthemodel

(>40–400kmdepth)wasconstructedﬁrst(seeFig.4UpperMantle

PartandFadeletal.(2013)fordetails).Adamsetal.(2012)used

a1kmthicklayeraboveandbelowthe40kmdepthtoconstrain

Mohodepthsandallowforlargevariationsinthevelocitybetween

thecrustalandtheuppermantlepartofthemodel.Hence,a1km

horizontallayeraboveandbelowadepthof40km(2kmthicklayer

betweendepths39–41km)wasconstructedtorepresenttheMoho

asintheoriginalmodel.Thisisinlinewithactualobservationson

theMohowhichoftenindicateitisalayerwithacertainthickness

andnotadiscreteboundary.Subsequently,thecrustalpartofthe

modelwasconstructedabove(Fig.4CrustalPart).

ThetomographicmodeldidnothavearealisticMoho

(4)

Fig.2.3DOOAresultsofthemodel.Thecrustalpart(0–40km)showstheOOAresultsoftheﬁveunits(U1–5).Theuppermantlepart(>40–400)showstheOOAresultsofthe differentuppermantleobjects:craton(C),shallowhighvelocity(SHV),rift(R),boundaryshallowhighvelocity(BSHV),boundaryshallowlowvelocity(BSLV),lowvelocity (LV),innerlowvelocity(ILV),boundarydeeplowvelocity(BDLV),anddeephighvelocity(DHV).Thereadershouldbeawarethattheobjectshaveareverseorientationin theEast-Westdirection(comparethemwithFig.4)whichisavisualizationdefectof3DOOAresults.

AdaptedfromFadeletal.(2015).

crustalandtheuppermantlepartsremainedthesame;however,

theﬁxedMohowasreplacedwithagravitybasedcrustal

thick-nessmapofAfrica(Tugumeetal.,2013).The2kmthickﬂatlayer

representing top and bottom of theMoho transitionzone was

modiﬁedtodescribetheundulationofthenewMoho.Nowa

real-isticMohotopographywasincluded,representingspatialdensity

variationsduetoMohoundulations.TheembeddedMohomodel

showedalargevariationincrustalthicknessesthatrangedfrom

∼28kmatcoastalregionsto∼42kminthecentralpartofthestudy

area(Fig.5).Thealternativeoptiontousecrustalthicknessfrom

CRUST2.0,acrustalthicknessmapthat wascreatedby

interpo-latingthecrustalvelocityestimatesfromactive-sourceproﬁling

techniquesandgeologicallyestimatedterrainages(Bassinetal.,

2000),wasnotconsideredbecauseofthedebateofusingCRUST2.0

indatapoorenvironments(seediscussioninvanderMeijdeetal.,

2013).

4.2. Object-basedinversion

Afterthedeﬁnitionofthetwoscenariomodelsbasedonthe

resultsofthe3DOOA,densityvaluesforeachobjectwereneeded

tocalculatethegravitysignalandcompareitwiththemeasured

satellitesignal.Schmidtetal.(2011)introducedalinear

object-basedinversionapproachinIGMAS+,basedontheminimummean

squareerror(MMSE)algorithmintroducedbyHaase(2008)based

ontheworkofSaether(1997).Theinteractiveobject-based

mod-ellingapproachaimstominimizethemisﬁtbetweenthecalculated

andthemeasuredgravitysignaluptoapproximationerror:

Calculated=measured−error=

n

i=1

iPi (1)

wherenisthetotalnumberofbodies,iistheconstantdensityof

bodyi,andPirepresentsthepolyhedrongravityeffectofbodyi(for

density1).

Eq.(1)simpliﬁesthesubsurfacestructureintoasmallnumber

ofobjects.Eachobjecthasaconstantdensityvalue(i)thatisthe

commonvalueofallvoxelswithinthisobjecti.Thecalculated

grav-ityeffectofobjectiwillbetheproductofitsdensityimultiplied

bythepolyhedrongravityeffectincaseofaunitdensityvalue(e.g.,

1t/m3_or₁_kg/m3_or_any_other_unit)._The_total_calculated_gravity

sig-nalwillbethesummationoftheresultingpolyhedroneffectsfor

thetotalnumberofobjectsn.Theobject-basedinversionapproach

(5)

Fig.3.3Dhistogramsofthemodel.(A)Thehistogramisthefull3Dhistogramofthemodelfrom0to400kmdepth.Theblackarrowat40kmdepthindicatesthesharp discontinuityinthemodelbecauseoftheMoho.(B)Thehistogramshowsthecrustalpart(0–40km)ofthemodel,whiletheblack-whitelinesindicatethethresholdvelocities forthedifferentunits.

insteadofthetraditionalinteractivetrialand errorapproaches.

Theestimatedoptimizeddensityvaluesi,directlyrelatedtothe

deﬁnedobjectsproduceacalculatedgravitysignalwithaminimum

misﬁt(error)tothemeasuredsignal.Itisalsopossibletoadd

con-straintstotheinversionprocessincaseofthepresenceofapriori

knowledgeaboutthedensityvaluesofsomeobjectsinthemodel.

Inthatcase, theobjectswithaprioriknowndensityvalues are

providedintheinversionalgorithmasconstantswhilethe

inver-sionalgorithmisonlyappliedtotheobjectswithunknowndensity

values,asshowninEqs.(2)and(3).

Calculated=

ip ipPip+

np npPnp (2)

HencebysubstitutionusingEq.(1),

Measured−error−

np npPnp=

ip ipPip (3)

whereindexipindicatestheobjectsthatwillbeinverted,andnp

indicatestheobjectsthatwillnot.Apointofattentionisthe

set-tingofthebackgrounddensityvalue.TheBougueranomalyisthe

resultofthedensitycontrastbetweenthebackgrounddensitiesof

theEarthandtheexistingdensitiesofthedifferentobjectsinthe

subsurface.IGMAS+allowsonedensitybackgroundvalueto

cal-culatedensityanomaliesand,hence,gravityanomaliesrelativeto

it.However,formodelswithlargeverticalandspatialdensity

con-trastsdownto400kmthebackgrounddensityvalueschangeswith

depth.Therefore,theuseofonevalueasadensitybackgroundis

notideal.Inordertodealwiththislimitation,relativeobject

den-sitycontrastswereusedwithrespecttoabackgrounddensityof

0.Theresultsfromtheinversionprocessarethusdensitycontrast

valuesbetweenobjectsinthemodelandthebackground

densi-tiesattheirequivalentdepths.Forexample,ifweassumedthatthe

densitybackgroundfor0–10kmis2.67t/m3_,_the_absolute_density

valueofunit1is2.453t/m3_(U1_lies_in_a_depth_range_0–10_km_with

adensitycontrastof−0.217,asshowninTable2).Themeasured

gravitysignalwasﬁlteredtoremoveshortwavelengthsinorderto

ﬁttheresolutionoftheinputdata(i.e.,theseismologicalmodel).

Thesmallestobjectinthemodelis approximately100km.This

requiredtheapplicationofalow-passcosineroll-offﬁlteronthe

measuredgravitysignal(Fig.6D).Thisremovedthehighfrequency

content,neededtobecomparablewiththecalculatedgravity

sig-nalbasedontheseismologicalmodel.Theﬁlteredgravitysignal

wasusedtoestimatetheoptimizeddensitycontrastvaluesofthe

differentobjects,andtoevaluatethedegreeofcorrelationwiththe

(6)

Fig.4. ThereconstructionprocessoftheextractedobjectsinIGMAS+.(FM)ThefullmodelimageshowstheimportedsectionsresultedfromOOAandthemeasuredgravity onthetopofthem.Theotherimagesshowthereconstructionprocessofeachobject:theﬁveunitsofthecrustalpart(U1–5),craton(C),shallowhighvelocity(SHV),rift(R), boundaryshallowhighvelocity(BSHV),boundaryshallowlowvelocity(BSLV),lowvelocity(LV),innerlowvelocity(ILV),boundarydeeplowvelocity(BDLV),anddeephigh velocity(DHV).

Table1showstheresultsoftheinversionprocessforthe

differ-entpartsofthetwomodels,withboththeﬂatandtheundulated

Moho,comparedtoEIGEN-6c.Inbothcases,thecalculatedgravity

signalfromtheuppermantlepartshowed69%correlationwiththe

measuredsatellitegravity.IncaseoftheﬂatMohomodel,theadded

Mohodiscontinuitylayerandthecrustalpartonlyhadaslighteffect

onimprovingtheresults(maximumof3%improvement).Incase

ofthecompletemodel,however,usingtherealisticMohosurface

basedonthecrustalthicknessmapofAfrica(Tugumeetal.,2013),

theresultsshowedamajorimprovementuptoa95%ﬁt.The

undu-latingMohodiscontinuityclearlyhasastronginﬂuenceontheﬁnal

ﬁt,whichistobeexpectedforalarge-scalestronganomalythatalso

hasasigniﬁcantcontributiontothegravitysignalamongall

mod-elledobjects.Anotherreasonfortheincreasedﬁtisthattheapplied

Mohotopographyisalsoderivedfromgravitydata.Therefore,it

isexpectedtogiveahighcorrelationandconsequentlycontribute

strongly.Thisstrongcontributionisreﬂectedinthe,unrealistically,

highrelativedensitycontrastvalue5.157t/m3_of_the_Moho_layer

thatresultsfromtheinversion,asshowninTable2.Thisisalso

partlybecausetheMohoismodelledasathinlayer(2kmthick)

andthatalltheuncertaintiesinthestrongnegativecrustal

den-sityvalues(U1−0.217t/m3_,_U2_−0.244_t/m3_,_U3_−0.278_t/m3_,_U4

−0.312t/m3_,_and_U5_−0.373_t/m3_,_Table₂₎_seem_to_accumulate_in

thehighpositivedensitycontrastvalueoftheMoho.

Table2(secondcolumn)showsthedensitycontrastvaluesfor

eachindividualobjectinthemodelresultingfromtheobject-based

inversionusingEIGEN-6c.Thevaluesforthecrustalpart(U1–U5)

ofthemodel,rifts(R),thelowvelocityzone (LV),and itsinner

part(ILV)showednegativedensitycontrast.Fortherifts(R),the

lowvelocityzone(LV),andtheinnerlowvelocityzone(ILV),the

valuesarematchingwiththelowvelocityvaluesfromthe

tomo-graphymodel.However,theinnerlowvelocity(ILV)objectshowed

Fig.5.ThecrustalthicknessmapofthestudyareaadaptedfromTugumeetal. (2013).

(7)

Fig.6.Thecalculatedgravitysignalbasedontheinversionprocess.(A)The100kmﬁlteredgravitysignal.(B)Thecalculatedsignalfromtheuppermantlepartofthemodel (>40–400km)with69%correlationwith(A).(C)Thecalculatedgravitysignalusingthefullmodelwith95%correlationwith(A,B,andC).(D)Thedifferencebetweenthe measured(A)andthecalculated(C)signal.(E)Thehistogramofdifferencemap(D).

lessnegativedensitycontrastthanthelowvelocity(LV)object,

whichcontradictstheseismictomographymodel,sincetheinner

lowvelocity(ILV)hadthelowestvelocity.Thecraton(C),the

shal-lowhighvelocity(SHV),andthedeephighvelocityzone(DHV)

objects,showedpositivedensitycontrastvaluesthatagreedwith

thehighvelocityinthetomographicmodel.However,thecraton

(C)hadahighervelocityinthetomographicmodelthanthe

shal-lowhighvelocity(SHV)object.Thiswascontrarytothedensity

Table1

Theinversionresultsforthedifferentpartsofthetwomodels.

Parameters Uppermantlepart UppermantlepartwithMoho Completemodel CompletemodelwithoutMoho FlatMohomodel

Correlation 0.69 0.71 0.72 0.70

Standarddeviation×10−5_m/s−2 ₂₅ ₂₄ ₁₄ ₂₅

UndulatedMohomodel

Correlation 0.69 0.71 0.95 0.89

(8)

Table2

Thedensitycontrast(t/m3₎_estimated_from_the_inversion_of_the_complete_model_with_the_undulated_Moho_using_different_gravity_models.

Object Eigen-6c(t/m3₎ _Eigen-6c2_(t/m3₎ _GOCO03S_(t/m3₎ _{GO-CONS-GCF-2-DIR-R4}_(t/m3₎ _{GO-CONS-GCF-2-TIM-R4}_(t/m3₎

Unit1(U1) −0.217 −0.215 −0.211 −0.211 −0.211 Unit2(U2) −0.244 −0.243 −0.24 −0.241 −0.241 Unit3(U3) −0.278 −0.277 −0.275 −0.275 −0.275 Unit4(U4) −0.312 −0.312 −0.308 −0.309 −0.308 Unit5(U5) −0.373 −0.374 −0.371 −0.372 −0.371 Craton(C) 0.008 0.008 0.007 0.007 0.007

ShallowHighVelocity(SHV) 0.011 0.011 0.01 0.01 0.009

BoundaryShallowHighVelocity(BSHV) 0.005 0.005 0.004 0.004 0.004

Rift(R) −0.002 −0.002 −0.003 −0.003 −0.003

BoundaryShallowLowVelocity(BSLV) −0.015 −0.015 −0.014 −0.014 −0.014

LowVelocityZone(LV) −0.052 −0.053 −0.053 −0.053 −0.053

InnerLowVelocity(ILV) −0.002 −0.002 −0.003 −0.003 −0.003

BoundaryDeepLowVelocity(BDLV) −0.067 −0.067 −0.062 −0.06 −0.058

DeepHighVelocity(DHV) 0.172 0.173 0.171 0.172 0.17

MOHO 5.157 5.154 5.132 5.137 5.139

Reference 0 0 0 0 0

Correlation 0.948 0.949 0.948 0.948 0.947

Standarddeviation 11.021 11 11.002 11.062 11.089

contrast,whichwaslowerforthecraton.Thiscontradictioncan beexplainedbythelargeuncertaintyinthevelocityvaluesthatis lessthan0.2km/satdepths<250km,andnotmorethan0.3km/s atgreaterdepths(Adamset al.,2012).Theboundariesbetween

ShallowHighVelocity(BSHV),ShallowLowVelocity(BSLV),and

DeepLowVelocity(BDLV)objectsshowedrealistictransition

val-ues.Theﬁnalmapsofthecalculatedgravitysignalfortheupper

mantlepartandthecompletemodelwiththeundulatedMohoare

showninFig.6.

4.3. Testingdifferentgravitymodels

Thepreviousmodelling wasdoneusing a combined

EIGEN-6cmodelconstructedfromsatellitedataandothergravitydata

sources.Thisraisedthequestionofhowmuchthedifferent

con-tent,resolutionandqualityofdifferentgravitymodelsinﬂuence

theﬁnalmodelparameters.Forthispurpose,differentcombined

andsatellitedataonlymodelswerecompared.EIGEN-6c,

EIGEN-6c2(GOCE,GRACE,LAGEOSandsurfacedata), GOCO03S(GOCE,

GRACE,CHAMP,SLR),GO-CONS-GCF-2-DIR-R4(GOCE,GRACE,and

LAGEOS),andGO-CONS-GCF-2-TIM-R4(GOCEonly)wereusedin

thistest(formoreinformationontheGOCErelatedmodelsseevan

derMeijdeetal.(2015)).Themodelsweredownloadedwiththe

sameresolutionandsubjectedtothesameprocessing(terrain

cor-rection,Bouguercorrectionand100kmcosineroll-offﬁlter).The

resultsforthedifferentgravitymodelsaresummarizedinTable2.

TheresultsshowedthatEIGEN-6c2,themostrecentand

high-estresolutionmodel,givesthebestresults.Thiswasexpectedsince

themodelhasbeenshowntohaveahigheraccuracythananyother

gravitymodel(Försteetal.,2013).However,thedifferencesinthe

models(Fig.7)werenotsigniﬁcant.Inaddition,thecorrelations

areallhighandcomparable,andstandarddeviationresultsvary

little(Table2).Whencomparedtothestandarddeviationthatis

obtainedforthebestmodel(seeTable1)onecanconcludethat

choicesmadein thecreationofthemodel areofmuch greater

consequencethanthevariabilitybetweenmodels(at100km

reso-lution).

5. Discussion

Itisdifﬁculttoquantifyhowmuch3DOOAhasimprovedthe

modellingofsubsurface objects.Anexpert modelermighthave

cometoasimilarresult,ormaybeevenbetter.Themainadvantage

isthattheboundariesarebasedonphysicalrulesandknowledge,

andthattheyareobjectiveandrepeatable.Anadditionaladvantage

wasthatthetopandbottomoftheobjectsweredetermined,

some-thingthattraditionallyhasbeendifﬁcultduetostrongsmearing

and/orsmoothingintheverticaldirection.

Inthisexampletheboundariesoftheobjectswereﬁxed,and

theinversiononlyinﬂuencedtheinternalphysicalpropertiesofthe

objects.Thisis,ofcourse,alimitationandcanbegivenmoredegrees

offreedominfutureattempts.Themethodwedevelopedisidealfor

incorporatingaprioriinformationintheidentiﬁcationofobjects.

Weonlyusedgeologicalknowledgeonthemaingeologicaland

tec-tonicfeaturesintheareabasedonAdamsetal.(2012),andcrustal

thicknessfromTugumeetal.(2013).Moreaprioriinformationcan

potentiallybeincorporated,suchasvelocity–densityrelationships,

pointobservationsfrom,forexample,receiverfunctionanalysis,

and petrochemicalknowledgeonthephysicalpropertiesof the

variousdeﬁnedobjects.Thiswayphysicalpropertiescanbe

fur-therconstrained.Inourexperimentsweobservedthatforsome

objectparametershadnon-realisticdensitycontrastvalues(e.g.,

theMoholayer).Byaddingrangeswithinwhichparametersmay

vary,theresultislikelytobemorerealisticintermsofretrieved

physicalproperties.

Becauseofthepoorlyresolvedcrustalpartoftheseismological

model (Adams et al., 2012), there is no strictly deﬁned

struc-turalintegritybetweenthecrustalandtheuppermantlepartof

themodel.ThesharpMohoboundarydeﬁnedinthetomographic

modelisnotphysicallyrealisticinallplaces,andreplacementby

anothercrustalthicknessmodelcreatesanuncertaintyinthe

phys-icalcontinuationofthephysicalproperties,therebyaddingmore

uncertaintyonstructuralintegrity.

The31%misﬁtvaluebetweenthecalculatedandmeasured

grav-ity,basedontheuppermantlepartofthemodel,canbeexplained

duetotheabsenceofthecrustalpartinthemodel.Thedramatic

improvementinthecorrelationandstandarddeviationafteradding

theMohoaredoubtfulwhenconsideringtheunrealisticallyhigh

density contrastof the Moho thathad a highinﬂuenceon the

measuredsignal.Therefore,inordertoprovethenecessitytoadd

arealisticMohosurface,theinversionwascarriedoutexcluding

theMoholayer(adensitycontrast=0)fromtheinversionusingnp

indexmentionedinEqs.(2)and(3).Theresultsoftheinversion

showedasigniﬁcantimprovementof89%inthecalculatedsignal

(Table1).Thisunderscorestheimportanceofknowledgeonthe

shapeoftheMohosurface.Inaddition,itreﬂectsthenecessityfor

densitycontrastconstraintsofsomeobjects,sinceaminorchange

inoneobjectinthemodelcausedsigniﬁcantchangesinthe

den-sitycontrastofallobjects.Table2showsthatanegligiblechange

(9)

Fig.7. MapsshowthedifferencebetweencombinedsatelliteandsurfacedatagravitymodelEIGEN-6c2thatshowedthebestﬁtintheinversionresults(Table2)andthe differentgravitymodels.(A)ThedifferencebetweenEIGEN-6c2andthecombinedsatelliteandsurfacedatagravitymodelEIGEN-6c.(B)ThedifferencebetweenEIGEN-6c2and satelliteonlygravitymodelGOCO03S.(C)ThedifferencebetweenEIGEN-6c2andsatelliteonlygravitymodel(GO-CONS-GCF-2-DIR-R4).(D)ThedifferencebetweenEIGEN-6c2 andsatelliteonlygravitymodel(GO-CONS-GCF-2-TIM-R4).

densitycontrastvalues.Anotherindicatorforthehighsensitivity

ofthedensitycontrastvaluesisthatthevaluesarenotstablewith

changesintheextentofthemeasuredgravitysignal.Sometests

weredonetoremovetheedgeeffectinthecalculatedgravity

sig-nalbyreducingthemeasuredgravityextentby4,8,and12%ofthe

totalextentofthestudyarea.Theresultsofthedensitycontrastfor

eachtrialwerealsodifferent.Moreover,theedgeeffectwasstill

presentintheforwardsignalandthemisﬁtbetweenthemeasured

andthecalculatedgravitysignal(Fig.6CandD)sincetheforward

calculationofIGMAS+isbasedonﬂatEarthassumptionsanddoes

nottakeintoaccountthesphericalshapeoftheEarth.

Theuncertainty in the location of theextracted objects are

relatedtotheuncertaintyoftheseismictomographicmodel.Adams

etal.(2012)statedthattheuncertaintiesformajorfeaturesinthe

modelarerangingfrom25to50kmfortheshallowzoneofthe

uppermantle partofthemodel (40–250km)and canbeupto

75–100kmfordepths>250km.Theuncertaintyofthecrustalpart

ofthemodel(0–40km)wasnotdiscussedinAdamsetal.(2012)

butwasassumedtobethesameasforthedepthrange40–250km.

Theuncertaintiesarereﬂectedintheresultsoftheobjectbased

inversion. The correlation factor and thestandard deviation of

theinversionresultsshowthatforthecompletemodel,including

theundulatedMohodensitycontrasts5.157t/m3_,_the_correlation

was0.95 and the standard deviation of the misﬁt wasaround

11×10−5m/s2₍_Fig.₆_D_and_E)._The_misﬁt_can_be_explained_by_the

inherentuncertainty intheseismicmodel thatalsoincludethe

objectboundaries.Theinversionalgorithmonlyallowedchange

inthedensitycontrasts.Thatresultedinsomeobjectsdensity

con-trasttobeunrealistic.However,ifthealgorithmwasallowedto

changetheboundarieslocationwithincertainlimitsthatarelinked

totheuncertaintyoftheinputdata,maybeboththelocationof

theobjects’boundaries and thedensity contrastvalues willbe

improved.Thiscanbedonethroughimplementing insightona

prioriknowledgeaboutthepossibledensitycontrastrangesofthe

objectsthatcanpreventaddingmorenon-uniquenesstothe

prob-lemduetotheavailableinﬁnitecombinationofboundarylocations

anddensitycontrastvalues.

6. Conclusion

Inthisresearchwedevelopedamethodtoextract3Dobjects

(10)

objectswerethenreconstructedinanobject-basedgravity

mod-elling environment (IGMAS+), and the density contrast value

foreachobject calculatedusinganobject-basedinversion

tech-nique.Thecalculatedgravitysignalfromtheuppermantle part

(>40–400km) of the model produced a 69% correlation with

observedgravity.Thelowﬁtismainlyduetotheabsenceofthe

Mohoundulation surface,and thecrustalpart.Byadding aﬂat

Moholayerwitha2kmthicknessandthecrustalobjects,basedon

thetomographicmodel,thecalculatedgravitysignalonlyimproved

by3%.However,byaddingtheMohoundulationsurfacebasedon

thecrustalthicknessmapofAfricabyTugumeetal.(2013)and,

theintra-crustalobjectsfromthetomographicmodelfromAdams

etal.(2012),theﬁtofthecalculatedgravitysignalwiththeobserved

gravitysigniﬁcantlyimprovedupto95%.However,theresultsfor

theMohoandthecrustalpartneedfurtherveriﬁcation,and

fur-therimprovementsareexpectedbasedonothergeophysicaldata

sources,suchasreceiverfunctions.

Theusageofthedifferentsatellitegravitymodels,either

satel-liteonlymodelsorsatelliteandterrestrialmodels,onlyhadaminor

effectontheinversionresults.Uncertaintyintheinversionislarger

thanthevariabilityduetotheuseofdifferentmodels(allat100km

resolution).

Integrating3DOOAina3Dobject-orientedgravitymodelling

environment canleadtoa signiﬁcantimprovement in

convert-ingthevoxel-basedgeophysicalmodelsintoobject-basedmodels

based on simple and objective knowledge-based classiﬁcation

rules. In addition, the object-based inversion approach can be

improvedbyusingthedifferentavailablesegmentationalgorithms,

which,forexample,allowobjectstogrowandshrinkautomatically,

andhenceallowsomeﬂexibilityforautomaticadjustmentofthe

boundarylocation,andthustooptimizeﬁttingresults.

Acknowledgments

WethankDr.AubreyaAdamsatWashingtonUniversityfor

mak-ingtheseismictomographymodelofthestudyareaavailablefor

thisstudy.WewouldliketothankProf.Dr.Hans-JürgenGötzeand

hisresearchteamatKielUniversityfor providingtheacademic

licensefortheIGMAS+software.

AppendixA. Supplementarydata

Supplementarydataassociatedwiththisarticlecanbefound,in

theonlineversion,athttp://dx.doi.org/10.1016/j.jag.2013.11.004.

References

Adams,A.,Nyblade,A.,Weeraratne,D.,2012.Uppermantleshearwavevelocity structurebeneaththeEastAfricanplateau:evidenceforadeep,plateauwide lowvelocityanomaly.GeophysicalJournalInternational189(1),123–142.

Bassin,C.,Laske,G.,Masters,G.,2000.Thecurrentlimitsofresolutionforsurface wavetomographyinNorthAmerica.EOSTransactionsAGU81,F897.

Bezada,M.J.,Zelt,C.A.,2011.Gravityinversionusingseismicallyderivedcrustal densitymodelsandgeneticalgorithms:anapplicationtotheCaribbean-South Americanplateboundary.GeophysicalJournalInternational185(2),577–592.

Birch,F.,1961.Thevelocityofcompressionalwavesinrocksto10kilobars:2.Journal ofGeophysicalResearch66(7),2199–2224.

Blaschke,T.,2010.Objectbasedimageanalysisforremotesensing.ISPRSJournalof PhotogrammetryandRemoteSensing65(1),2–16.

Brocher,T.,2005.Empiricalrelationsbetweenelasticwavespeedsanddensity intheEarth’scrust.BulletinoftheSeismologicalSocietyofAmerica95(6), 2081–2092.

Cammarano,F.,Tackley,P.,Boschi,L.,2011.Seismic,petrologicalandgeodynamical constraintsonthermalandcompositionalstructureoftheuppermantle:global thermochemicalmodels.GeophysicalJournalInternational187(3),1301–1318.

Christensen,N.,Mooney,W.,1995.Seismicvelocitystructureandcompositionofthe continentalcrust;aglobalview.JournalofGeophysicalResearch-SolidEarth100 (B6),9761–9788.

Ebbing,J.,Braitenberg,C.,Götze,H.,2001.Forwardandinversemodellingofgravity revealinginsightintocrustalstructuresoftheEasternAlps.Tectonophysics337 (3–4),191–208.

Fadel,I.,Kerle,N.,Meijde,M.,2015.3Dobjectorientedimageanalysisofgeophysical data.GeophysicalJournalInternational35,44–53.

Förste,C.,Bruinsma,S.,Flechtner,F.,Marty,J.-C.,Dahle,C.,Abrykosov,O.,Lemoine, J.-M.,Neumayer,H.,Barthelmes,F.,Biancale,R.,etal.,2013.Eigen-6c2-anew combinedglobalgravityﬁeldmodelincludingGOCOdatauptodegreeand order1949ofGFZ-PotsdamandGRGS-Toulouse.In:EGUGeneralAssembly ConferenceAbstracts,vol.15,p.4077.

Förste,C.,Bruinsma,S.,Shako,R.,Marty,J.,Flechtner,F.,Abrikosov,O.,Dahle,C., Lemoine,J.,Neumayer,K.,Biancale,R.,etal.,2011.EIGEN-6-anewcombined globalgravityﬁeldmodelincludinggocedatafromthecollaborationof GFZ-PotsdamandGRGS-Toulouse.In:GeophysicalResearchAbstracts,vol.13.

Gardner,G.,Gardner,L.,Gregory,A.,1974.Formationvelocityanddensity:the diag-nosticbasicsforstratigraphictraps.Geophysics39(6),770–780.

Godfrey,N.,Beaudoin,B.,Klemperer,S.,1997.OphioliticbasementtotheGreat Val-leyforearcbasin,California,fromseismicandgravitydata:implicationsfor crustalgrowthattheNorthAmericancontinentalmargin.GeologicalSociety ofAmericaBulletin109(12),1536–1562.

Götze,H.,Lahmeyer,B.,1988.Applicationofthree-dimensionalinteractivemodeling ingravityandmagnetics.Geophysics53(8),1096–1108.

Haase,C.,2008.Improvedestimation ofsubsurfacemagneticproperties using minimummean-square error methods.Christian-Albrechts-UniversityKiel, Germany(Diplomathesis).

Korenaga,J.,Holbrook,W.S.,Detrick,R.S.,Kelemen,P.B.,2001.Gravityanomalies andcrustalstructureatthesoutheastGreenlandmargin.JournalofGeophysical Research:SolidEarth106(B5),8853–8870.

Liu,Q.,Gu,Y.,2012.Seismicimaging:fromclassicaltoadjointtomography. Tectono-physics566–567,31–66.

Martha,T.,Kerle,N.,vanWesten,C.,Jetten,V.,Kumar,K.,2011.Segment opti-mizationanddata-driventhresholdingforknowledge-basedlandslidedetection byobject-basedimageanalysis.IEEETransactionsonGeoscienceandRemote Sensing49(12),4928–4943.

Mishra,D.,Kumar,M.R.,Arora,K.,2012.Longwavelengthsatellitegravityandgeoid anomaliesoverHimalaya,andTibet:lithosphericstructuresand seismotecton-icsofdeepfocusearthquakesofHinduKushPamirandBurmesearc.Journalof AsianEarthSciences48,93–110.

Ravat,D.,Lu,Z.,Braile,L.,1999.Velocity–densityrelationshipsandmodelingthe lithosphericdensityvariationsoftheKenyaRift.Tectonophysics302(3–4), 225–240.

Saether,B.,1997.Improvedestimationofsubsurfacemagneticpropertiesusing min-imummean-squareerrormethods.Norgesteknisk-naturvitenskapeligeuniv, Trondheim,Norway(Ph.D.Thesis).

Sanchez,J.,Götze,H.,Schmitz,M.,2011.A3-Dlithosphericmodelofthe Caribbean-SouthAmericanplateboundary.InternationalJournalofEarthSciences100(7), 1697–1712.

Schmidt,S.,Götze,H.,Fichler,Ch.,A.M.,2010.IGMAS+:anew3Dgravity,FTGand magneticmodellingsoftware.Extendedabstract.Geoinformatik,57–63.

Schmidt,S.,Plonka,C.,Götze,H.,Lahmeyer,B.,2011.Hybridmodellingof grav-ity,gravitygradientsandmagneticﬁelds.GeophysicalProspecting59(6,SI), 1046–1051.

Schoenmeyer,R.,Prvulovic,D.,Rotarska-Jagiela,A.,Haenschel,C.,Linden,D.E.J., 2006.Automatedsegmentationoflateralventriclesfromhumanandprimate magneticresonanceimagesusingcognitionnetworktechnology.Magnetic Res-onanceImaging24(10),1377–1387.

Tugume,F.,Nyblade,A.,Juli,J.,vanderMeijde,M.,2013.Precambriancrustal struc-tureinAfricaandArabia:evidencelackingforsecularvariation.Tectonophysics,

http://dx.doi.org/10.1016/j.tecto.2013.04.027(inpress).

vanderMeijde,M.,Juli,J.,Assumpo,M.,2013.GravityderivedMohoforSouth Amer-ica.Tectonophysics,http://dx.doi.org/10.1016/j.tecto.2013.03.023(inpress). vander Meijde,M., Pail,R.,Bingham, R.,Floberghagen, R.,2015. GOCEdata,

models,andapplications;areview.Int.J.Appl.EarthObs.Geoinform. 35, 4–15.