Contribu ții la descifrarea structurii profunde din zona seismogen ă Vrancea și regiunile adiacente utiliz ând date gravimetrice și geodezice… [602865]

1
Universitatea București
Facultatea de Geologie și Geofizică

PROIECT DE CERCETARE

Contribu ții la descifrarea structurii profunde din zona
seismogen ă Vrancea și regiunile adiacente utiliz ând date
gravimetrice și geodezice

Coordonato r științific Doctorand: [anonimizat]. Marian Ivan Ing. Natalia -Silvia Chira
(Asimopolos )

– 2012 –

2

ANOMALIA FREE AIR PE TERITORIUL ROMÂNIEI

1) Anomalia free -air calculata pe baza valorilor medii de gravitate Bouguer

Deducerea anomaliei free air („ în aer liber”) folosind harta altitudinilor medii și
harta an omaliei Bouguer a suscitat interesul încă din anii 1955 , prin prezentările ținute de
Mircea Socolescu și D. Bișir la Institutului de Fizi că, care au fost publicate în volumele
de Comunicări ale A cademiei Române.
Acest lucru a afost de un real ajutor in pr ogresul geofizicii, datorit ă faptului ca în
țara noastr ă valorile de gravitate ob ținute în numeroasele lucr ări de teren se reprezentau
numai sub forma anomaliei Bouguer, ceea ce însemna un volum imens de lucru în
recalcularea valorilor de gravitate m ăsurat e pentru fiecare sta ție.
Obținerea hărții „în aer liber” prin procedeul grafic ( destul de simplu dar
neriguros definit ), se face pe baza reducerii hărții Bouguer cu ajutorul hărții altitudini lor
medii și are avantajul că ține seama și elimină influența ma selor laterale, echivalând,
oarecum , cu transformarea reliefului accidentat într -un relief ondulat.
Harta anomaliei „în aer liber” astfel obținută s -a putut utiliza în calculele
geodezice.

În determinarea unei solu ții gravimetrice pentru geoid, harta anomali ei free air a
fost obținută de către Ioane et al . în diverse etape ale unor studii pentru structura adanc ă
(1994,1995,1996) studii aflate în arhiva Fondului Geologic al IGR.
Astfel (Ioane 1994) au fost utilizate valori de gravitate în blocuri de 5’ x7,5’,
pornindu -se de la Harta anomaliei Bouguer sc 1 : 200 000, cu izolinii trasate la nivelul de
1 mgal, și apoi evalu ând direct valorile de gravitate pe foile de harta, în aceia și rețea, la

3
scara 1 : 100 000 , neținându-se cont de valorile punctuale cu distribu ție neordonat ă și cu
riscul calit ății hărții si al abilit ății opera torului . Totodat ă s-au utilizat ș i valorile medii ale
cotelor topografice în aceia și rețea.
Pentru calculul anomaliei free air s-a folosit formula dat ă de Frosberg în 1994
(Terrain effect s în geoid computation) , iar pentru calculul gravita ții normale s -a folosit
formula Cassinis.

Unul din impactele reducerii locale ale terenului este s ă se scadă corelația
anomaliei free air cu înălțimea, evitând astfel erorile ce apar cand punctele de măsură ale
gravității au tendinț a de localizare în munți și văi.

Harta gravității României a fost calculată la Institutul Geologic al Rom âniei
bazându -se pe 78000 puncte de observare (media de acoperire de 0,33 puncte pe km2) de
către Ioane, D., și R adu, I., 1995 (Global geopotential models and gravity data for the
territory of Romania . In: Su nkel, H., Marson, I. (Eds.), Gravity and Geoid. Joint IGC/ICG
Symposium Graz 1994, IAG Symposium No. 113. Springer, Heidelberg, pp. 640 –646,
http://adsabs.harvard.edu/abs/1995grge.conf..640I , lucrare citată în Blanka Sperner et al.
2004 ).
Formul ele de calcul a anomaliei free air din lucra rea amintita mai sus, ca și î n
lucrarea Ioane et Atanasiu, 1998 nu a putut fi consultate în studiul de față, ele nefiind la
Fondul Geologic al IGR .

Am pr eluat din Sperner et al, 200 4, harta anomaliei fr ee air ș i profilul NW -SE
din lucra rea publicat ă în Tectono phys ics, așa cum apar mai jos.

4

Fig. 1. Harta anomaliei Free -Air a Rom âniei (Ioane and Atanasiu, 1998). Din 78,000
puncte de observație, valoarea medie a fost calculată într -un grid dreptunghiular de 5 V
latitudine x 7.5 V longitudine. Steluța albă reprezintă zona seismogenă Vrancea (după B.
Sperner et al. 2004 ).

5

Fig. 2. Harta anomaliei free air Io ane 1998, digitizat ă pe izolinii în Surfer 08 (4000 de
puncte)

2) Anomalia free -air pe baza modelului de geopotențial EGM08
(http://bgi.omp.obs -mip.fr/ )

Formatul datelor gravimetrice și calculul anomaliilor conform national
geospat ial-inteligence agency

Gravitatea observat ă (sau măsurată) g este valoarea gravității din locația stației.
Toate valorile trebuie ajustate la International Gravit y Standardi zation Net din 1971.
Gravitatea teoretică (normală) γ este valoarea de referință a gravității obținută din
câmpul gravific al Wor ld Geodetic System (WGS84) elipsoidul de revoluție de
referință. Aceasta este dată de formula:

6

 
 2
21 0.0019318526241 sin
978032,53359
1 0.0066943799014 sinmGal




unde
 este latitudinea geodezică. In formă analitică această ecu ație este dată de:

2/12 22
sin 1sin 1

ek
e
, unde:
1
eP
abk
a = semiaxa mare (elipsoidul WGS84),
b = semiaxa mică (elipsoidul WGS84),
γP = gravitatea normală la poli (WGS84 EGM96 – Earth Gravity Model),
γe = gravitatea normală la ecuator (WGS 84 EGM96 – Earth Gravity Model),

= latitudinea geodezică,
e2 = pătratul primei excentricități (elipsoidul WGS84)

Corecția gravității atmosferice (δg A) este corecția care se adaugă gravității
observate. Acest lucru este necesar deoar ece constanta gravitațională pentru WGS84
include masa atmosferei. Ea este dată de:
1.047
0.11610000.87 0
0.87 0h
A
Ag e mGal pentru h
g mGal pentru h
  


unde h este înălțimea față de nivelul mării.

Gradientul vertical al gravității normale (
) este rata de schimbare a gravită ții
teoretic e pe direcția verticală la supr afața elipsoidului. El este dat de relația:
 2sin 2 1 2 f mfa h

Gradientul vertical de ordinul 2 poate fi apreciabil la înălțimi h mai mari față de
nivelul mării:
a h622


7
Anomalia Free -Air (Δg f) este defin ită ca diferența dintre gravitatea observată pe
suprafața fizică (P) și gravitatea normală pe teluroid (Q). Teluroidul este definit ca
suprafața unde potențialul gravității normale este egal cu potențialul actual pe
suprafața fizică. Înălțimea deasupra eli psoidului la care potențialul normal este egal
cu potențialul real pe suprafața fizică este numită înălțime normală. Formulele
anomaliei gravității din figurile urm ătoare au la bază considerentul că în ălțimea
normală este egală cu altitudinea stației de mă sură a gravității.

Calculul anomaliilor Free -Air ( Heiskanen & Moritz – Physical Geodesy, 1967, p.
293) rezultă din: Δgf =gP – γQ ,
unde gP este gravitatea reală măsurată pe suprafața fizică iar γQ este gravitatea
normală pe suprafața teluroidului:

2*
22
*
0!21HhHhQ 




 
,
unde:
0 este gravitatea normală pe ellipsoid iar
*H este înălțimea normală.

Simboluri si notații:

Δgf anomalia gravității Free -Air (mgal)
δgA corecția atmosferică (mgal)
γ gravitatea teoretică (mgal)
g gravit atea observată (mgal)
h elevația punctului de observație (metru)
(suprafața terenului, a apei, a gheții)
d elevația suplimentară (adâncimea oceanului, (metru)
lacului, gheții sau a instrumentului
ρ densitatea (gram/cm3)

8

Paramet ri elipsoidului WGS 84 folosiți în ecuații sunt:
Semiaxa mare a = 6 378 137 m
Semiaxa mică b = 6 356 752.3142 m
Excentricitatea e = 0.081819190842622 , e2 = 0.00669437999014
Viteza unghiulară ω = 7 292 115 * 10-11 radians / sec ± 0.1500 * 10-11 radians / sec
Aplatizarea f = 0.00335281066474
Gravitatea ecuatorială normală γ e = 9.7803253359 m / sec2
Gravitatea normală la poli γ P = 9.8321849378 m / sec2
Constanta gravității normale : k = 0 .00193185265241
0684 0034497865.022
GMbam

Constanta gravitațională
2 3 8 2 3 8sec/ 101.0 sec/ 10 418.004 9863 m m GM   

Harta anomaliei free air a fost obținută prin interpolarea gridul ui datel or furnizat e
de la BGI (BUREAU GRAVIMETRIQUE INTERNATIONAL ), cu un numar de 28512
de valori ale anomaliei notate EGM08_FAA_Global_2190 _2,5 m . Harta anomaliei free
air pe teritoriul României are la bază și măsurători terestre (123 puncte).

9

Fig. 3. Date gravimetrice terestre în evidența BGI

Pe hart a anomaliei free air am trasat profilul NW -SE exact ca în harta anomaliei
free air descr isa la punctul 1 (Ioane, 1998) , apoi am construit graficul anomaliei FAA pe
acest profil.
Modelul geopoten țial EGM08 (Earth Gravitational Model 2008) a fost lansat
public de către Național Geospațial -Intelligence Agency (NGA) (N. K. Pavlis, S. A.
Holmes, S . C. Kenyon and J. K. Factor. “An Earth Gravitational Model to Degree 2160:
EGM2008”. EGU General Assembly 2008, Vienna, Austria, April 13 -18, 2008) , se
dezvoltă până la gradul 21 90, si conține coeficienți de extindere la 2159 al armonicelor
sferice (compl etează vechiul model EGM96) și reprezintă un progres, dorindu -se a se
realiza o bază de date complete cu valori ale anomaliei gravității .
Dezvoltarea EGM 2008 utilizează bazele de date la nivel mondial cu cele mai
bune informaț ii disponibile despre gravit atea terestră , marină , precum ș i sursele GPS și
de la misiunile prin satelit (Grace) , și care să conțină și cele mai recente date și metode de
modelare pentru sisteme maritime și terestre. Pentru ondulațiile geoidului acest model se
referă la elipsoidul de referință WGS84.

10

Fig. 4. Harta anomaliei Free Air (după date de la BGI) interpolate cu metoda kriging și
netezite cu funcții spline. Matricea de grid are 466 linii și 598 coloane.

11

3) Comparație între cele doua grafice ale anomaliei FAA pe profilul NW-SE,
din cele doua hărți ale anomaliei free air (BGI și Io ane 1998)

Fig. 5 Diferența dintre valorile de pe profilul reprezentat pe cele două hărți Free Air .

-80-60-40-200204060
22.952572 23.151181 23.326873 23.433816 23.578954 23.808117 24.182419 24.434499 24.785885 25.144909 25.374073 25.587959 25.809484 25.946982 26.260173 26.367117 26.451144 26.588642 26.848361 27.085164 27.390717 28.001821 28.887922LongitudineFree Air (mgal)
0.880.90.920.940.960.981
BGI
Ioane
CorelatieCoeficient
de
corelatie

Fig. 6. Reprezentarea celor două profile și corelația dintre ele . Pe ordonata din stânga
sunt reprezentate valorile de pe cele două profile (în mgal) iar pe ordonata din dreapta
este reprezentată variația coeficientului de corelație, folosind o fereastră mobilă d e 9
eșantioane și pasul de 1 eșantion.

Se poate observa că acest coeficient de corelație are valori cuprinse între 0. 89 și 0,999, iar
în jumătatea de profil din sud¤estul hărții acest coeficient are valori de 0,99, ceea ce
reprezintă o corelație aprope p erfectă.

12

Fig. 7 . Harta anomaliei free air BGI (roșu) suprapusă peste harta anoma liei free air Ioane
1998 (negru).

13

Fig. 8 Harta diferențelor pentru anomalia Free Air între BGI și Iane și Atanasiu(1998).

Am realizat gri dul într -o matrice de 5000 puncte (50×100) folosind metoda Kriging și
nivelarea cu funcții spline, pentru ambele hărți.
Am făcut apoi diferența dintre cele două griduri (Fig. 8)
Am calculat indicele de corelație dintre cele două hărți pentru toate cele 500 0 de puncte
din gridurile lor, obținând valoarea 0,886, ceea ce reprezintă un indice bun de corelație.
Dacă iau în considerație rezoluția grafică a hărții Free Air din (Sperner, 2004), pe care am
scanat -o am constatat că în zonele cu gradient mare și impli cit cu izolinii foarte dese,
scanarea nu a putut fi foarte precisă, datorită neclarității grafice a izoliniilor.
După eliminarea din hartă a zonelor cu gradient foarte mare (ce mulează linia carpaților),
am recalculat indicele de corelație, obținând 0,99 c eea ce reprezință un indice de corelație
foarte bun.
Deci se poate afirma că harta anomaliei Free Air obținută din datele de la BGI respectă la
fel de bine toate particulari tățile ca și harta realizată din datele la sol (Sperner, 2004).

Prin eliminea a 10 % din puncte apropierea dintre cele două hărți crește considerabil:

14

Fig. 9 Harta diferențelor pentru anomalia Free Air între BGI și Iane și Atanasiu(1998),
din care au fost eliminate 10% din punte. Izoliniile mai groase sunt izoliniile de valoare
zero.

Distribuția diferențelor dintre cele două hărți într -o histogramă are forma clopotului lui
Gauss, în care peste 50 din măsurători se încadrează în intervalul cuprins între -1 mgal și
+1 mgal.

Coeficien tul de corela ție a celor două seturi de date (Free Air de la BGI și Free
Air din Ioane, 1998) folosind factorul de corelație pentru o serie bidimensională de date:
  2 2 2 2/) ( ) () )( (

 
i i i ii i ii
xyy yn x xny x yxnr

unde: x i, yi sunt cele două serii de date; i = {1,…,n}, n = numărul de date din
aceaste două serii.
Factorul de co relație l -am calculat pentru o fereastră mobilă ce conține 9 perechi
de valori (x i, yi). Această fereastră mobilă am deplasea t-o, cu pasul de un eșantion, pe

15
întregul semnal. Astfel, factorul de corelație l-am calculat în 60 de puncte din cele 68 ale
profi lului; am eliminat primele 4 și ultimele 4 puncte din cele 68 de puncte ale profilului.

-80-60-40-200204060
22.952572 23.151181 23.326873 23.433816 23.578954 23.808117 24.182419 24.434499 24.785885 25.144909 25.374073 25.587959 25.809484 25.946982 26.260173 26.367117 26.451144 26.588642 26.848361 27.085164 27.390717 28.001821 28.887922LongitudineFree Air (mgal)
0.880.90.920.940.960.981
BGI
Ioane
CorelatieCoeficient
de
corelatie

Factorul de corelație pentru întraga serie de 0.995639, iar variația coeficientului de
corelație cu fereastră mobilă cu serii formate din câte 9 valori s -a încadra t în intervalul
0.89 și 1

In studiul de față am testat două metode de interpolare Kriging și metoda minimei
curburi, urmată de o netezire cu metoda funcțiilor spline.
Metoda Kriging, dezvoltată de D.G.Krige în 1950 pentru aplicații în geologie, este
meto da cel mai des folosită pentru interpolare pentru că se bazează pe formularea
statistică a celei mai bune estimații liniare. Un concept important pentru obținerea acestei
metode este bazat pe variograma experimentală:
   hjYjXiYiXji hN
hNjZiZhNCh      , , , )( ,
)() ())((1
21
,
C(N(h)) este numărul elementelor setului N(h). Variograma empirică este aproximată
folosind modele teoretice cum ar fi:

16
Modelul liniar
0 , 0 )(  h Sh Ch
Modelul gaussian







22
exp1 0 0 )(
ahCC Ch
Modelul exponențial




ahCC Ch exp1 0 0 )(
Metoda Kriging este folosit ă pentru estimarea interpolării funcției f(x,y) în punctul (X,Y)
creată special pentru aplicații în geologie, conservând bine trendul în extrapolări dar
având dezavantajul că poate produce izolinii circulare în jurul câte unui punct.
Metoda minimei curburi a fost dezvoltată de WHF Smith și P.Wessel (a se vedea [1]), în
1990 . Suprafața interpolată generată prin această metodă este similară cu o trecere a unei
lamele subțiri elastice prin fiecare din valorile date astfel încât se obține o sumă de
minimă îndoi re.Algoritmul metodei minimei curburi este bazat pe soluția numerică a
ecuației diferențiale biharmonice modificate:
0 ,2,41    yxf Tyxf T
cu trei condiții la limită:
 022
1 
nfT
nfT
,
02

nf pe muchii ,
02

yxf pe colțuri,
unde: T este parametrul de tensionare și are valori cuprinse între 0 și 1,
n este normala la limită,
22
222
y x
 este operatorul Laplacian iar
224
244
44 22 4
yx y x 


este operatorul biarmonic.
Această metodă este foarte bună pentru câm purile potențiale dar are dezavantajul că nu
extrapolează realist trendul izoliniilor în zonele fără înregistrări.
În cazul de față , după testarea acestor două metode, am exemplificat metoda de
interpolare Kriging, urmată de o netezire cu metoda funcțiilor spline .

Bibliografie

17
– Blanka Sperner, Dumitru Ioane, Robert J. Lillie – Slab behaviour and its surface
expression: new insights from gravity modelling in the SE -Carpathians
Tectonophysics 382 (2004) 51 – 84

– Miroslav Dressler – Art of Surface Interpola tion – PhD Thesis Kunštat, 2009

– Socolescu M, Bisir,D – Harta anomaliilor izostatice și ,,în aer liber ;pentru
Oltenia și Muntenia, Revista de studii si cercetari de fizică, vol VII, nr 2, 1956

– Biroul Gravimetric Internațional (http://bgi.omp.obs -mip.fr/ )