Budapesti Mőszaki és Gazdaságtudományi Egyetem, Építımérnöki Kar 2005. évi Tudományos Diákköri Konferencia www.vit.bme.hu/tdk/2005
Geopotenciális modell számítása GRACE bázisvonal alapján Szerzı: Paizs Zoltán,
[email protected], V. évf. építımérnök hallgató Konzulens: Dr. Földváry Lóránt, MTA-BME tudományos munkatárs, BME Általános- és Felsıgeodézia Tanszék Összegzés A GRACE egy gravitációs tér meghatározására irányuló NASA/GFZ projekt, amely 2002 tavaszától várhatóan 5 éven keresztül szolgáltatja a gravitációs mezı idıbeli és térbeli változásairól információt. A dolgozat keretében az energia integrál segítségével 4 hónap hosszúságú GRACE észlelésekbıl egy geopotenciális modellt számoltunk. A feldolgozáshoz a GRACE mőholdak ún. „1b szintő” (Level 1b) adatait használtuk. A számítások során mindvégig kinematikus pályákat használtunk fel. A megoldás a két mőhold közé létesített mesterséges bázisvonalra felírt energia-megmaradás törvény alapján történt. A pályaadatokból a bázisvonal mentén potenciálváltozást határoztuk meg, és ez alapján egyenlítettük ki a földi potenciál gömbfüggvénysorának együtthatóit. A bázisvonal alapján végzett meghatározás nem használja fel a mőholdak közötti, nagypontosságú mikrohullámú méréseket (KBR), így pusztán geometriai adatokból, a mőholdak pályája alapján kaptunk egy a mikrohullámú méréseket felhasználó megoldással analóg eljárást. Ez adatfeldolgozási szempontból egy jó kiindulópontja a KBR méréseket is felhasználó megoldásnak, melyet a jövıben szeretnénk kidolgozni. A dolgozat egy hosszabb távú, nagyobb szabású tanszéki projekt része.
Gravity model computation based on GRACE baseline Abstract GRACE is a NASA/GFZ joint gravity satellite project, which provides temporal and spatial gravity variations for 5 years from spring of 2002. In the frame of this study 4 months of GRACE data has been processed using the energy integral, and a low-degree gravity model has also been derived. Level 1b data has been employed, for gravity inversion kinematic orbit was used. The energy integral was defined on a GRACE inter-satellite baseline. The observation vector for the adjustment of spherical harmonic coefficients contains potential differences of the satellite. This solution makes no use of the KBR measurements, though it provides a nearly identical formulation to that. Therefore it is a good point of starting. The study is part of a larger project at the Department of Geodesy and Surveying at BUTE. Budapest, 2005. október 26.
Tartalomjegyzék 1
Alapok ........................................................................................................................... 3 1.1 Mesterséges holdak mozgása ................................................................................ 3 1.2 A GRACE mőholdak............................................................................................. 4 1.3 Perturbációk........................................................................................................... 4 2 A közvetítıegyenlet levezetése ..................................................................................... 6 2.1 Energia-integrál ..................................................................................................... 6 2.2 Low-Low Satellite to Satellite Tracking egyenlete ............................................. 11 2.2.1 Low-low SST bázisvonalon alapuló megoldásának közvetítıegyenlete..... 11 2.2.2 Low-low SST mikrohullámú távolságmérésen alapuló megoldásának közvetítıegyenlete ....................................................................................................... 12 3 Mérések feldolgozása .................................................................................................. 15 3.1 Adatok ................................................................................................................. 15 3.2 Pályatípusok ........................................................................................................ 15 3.2.1 Kinematikus pálya ....................................................................................... 15 3.2.2 Dinamikus pálya .......................................................................................... 15 3.2.3 Féldinamikus pálya (reduce dynamic orbit) ................................................ 16 3.3 Programok ........................................................................................................... 16 3.3.1 Sebességszámítás......................................................................................... 17 3.3.2 Adatbeolvasás és feldolgozás ...................................................................... 17 3.3.3 Energia-integrál kiszámítása........................................................................ 20 3.3.4 GRACE-EGM potenciálkülönbség meghatározása a drift levonásához ..... 23 3.3.5 A hosszú hullámhosszú tag (drift) meghatározása és levonása................... 25 3.3.6 A kiegyenlítés és a kapott eredmények elemzése........................................ 26 3.4 Eredmények ......................................................................................................... 29 4 Mellékletek .................................................................................................................. 30 4.1 Rövidítések .......................................................................................................... 30 4.2 Irodalomjegyzék .................................................................................................. 30 4.3 Feldolgozott adatok formátuma........................................................................... 31 4.4 A mőholdak pillanatnyi helyzete......................................................................... 33
2
1 Alapok 1.1 Mesterséges holdak mozgása A bolygók mozgásának törvényszerőségeit már Kepler felismerte és tapasztalati úton felállított törvényeivel leírta. Newton az általános tömegvonzás törvényének meghatározása után, ennek alapján megadta a Kepler-törvények fizikai magyarázatát is, mely szerint bármely égitest (Föld) közelében mozgó égitestek (mesterséges holdak) mozgását az égitest nehézségi erıtere szabja meg. Tehát ha megfigyeljük a mesterséges holdak mozgását, akkor lehetıségünk nyílik arra, hogy ezzel meghatározzuk, esetünkben a Föld nehézségi erıterét a GRACE A és GRACE B mőholdak mérései alapján. A mesterséges holdak mozgásának vizsgálata azonban sok szempontból eltér egy szabályos Kepler-féle mozgásától. Ennek oka elsısorban az, hogy a mesterséges holdak pályaszámításához sokféle perturbáló erıhatást kell figyelembe venni. Míg számos esetben az égitestek tömegpontoknak tekinthetık, addig a mesterséges holdak mozgása során – ezek viszonylagos felszíni közelsége miatt – a gravitációs perturbációk kiszámítása során az égitestek nem gömbszimmetrikus tömegeloszlásból származó hatásokat is figyelembe kell venni. Emellett a mesterséges holdaknál a gravitációs perturbációk mellett jelentıs hatással bírnak a nem-gravitációs effektusokat is. Mint minden földi mesterséges holdnál a fı gravitációs perturbációk magától a Földtıl származnak, melynek legnagyobb tagja a Föld egyenlítıi lapultsága. A nemgravitációs perturbációk közül, pedig a légköri fékezésbıl származó perturbáció a leginkább számottevı az ilyen alacsony magasságban keringı mőholdak esetén. Emellett a sugárnyomás perturbációkat sem szabad elhanyagolni esetünkben. A mőholdak mozgását perturbáló erık befolyásolják, melynek fı tagja természetesen a Föld, mint tömegpont által kifejtett tömegvonzási potenciál. Csak ezt az erıt figyelembe véve a mozgásegyenlet a következıképpen alakul:
&r& = −
kM r r2 r
A Föld valóságos gravitációs erıtere nem centrális erıtér. Ennek eltérését a potenciál gömbfüggvény-sorával írhatjuk le, melynek eredményeként a mozgásegyenlet a következıképpen fog alakulni:
&r& = −
kM r + grad R(J n , Clm , Slm , r) r2 r
A földi nehézségi erıtér leírása hagyományos geocentrikus, Földhöz kötött koordinátarendszerben történik. Ebben a rendszerben a Föld forgása keltette centrifugális erı is állandóan jelen van: &r& = −
kM r + grad R(J n , Clm , S lm , r) + a c r2 r
A további kisebb erıhatásokkal az 1.3-as pontban foglalkozunk. Az egyenlet megoldása során most már olyan Kepler-féle pályaelemeket kapunk, amelyek a Föld valódi nehézségi erıterét írják le az idı függvényében. A Föld valóságos nehézségi erıterében mozgó mőholdak (tömegpontoknak tekinthetık a Föld méreteihez képest) pályája bonyolult térgörbe lesz, melyek elemi szakaszai idıben változó mérető, alakú és helyzető Kepler-féle ellipszis pályák elemi darabjainak tekinthetık.
3
1.2 A GRACE mőholdak Napjainkban a gravitációs tér mérésénél egyre inkább elıtérbe került az őrgravimetria. Elsıként a német GFZ CHAMP mőholdját állították pályára 2000-ben, majd ezt követıen 2002-ben a NASA és a GFZ közös vállalkozásával a GRACE projekt keretében két mőholdat lıttek fel. Illetıleg az ESA gradiometriai mőholdja (GOCE) várhatóan 2006-ban indul útjára. A GRACE projekt két mőholdat takar, melyek Föld körüli pályán mintegy 485 km magasságban a Föld felszíne felett, közel poláris pályán (melynek névleges hajlása 89°), egymást 220 ± 50 km távolságban követve keringenek. A GRACE mőhold pár rendelkezik GPS vevıkkel, melyek segítségével folyamatosan mérik a pozíciójukat. Ez egy magas-alacsony mőholdkövetésnek (high-low SST) felel meg. A mőholdak pályaadataiból lehetıség nyílik arra, hogy egy bázisvonalat állítsunk fel. A pozíció meghatározásának a pontossága itt ±3 cm alatt van, tehát a két GRACE mőhold közötti távolság-meghatározás pontossága ±5 cm alatt van. Természetesen a mőholdok folyamatosan mérik a közöttük lévı távolságot két frekvencián (K 24GHz és Ka 32 GHz-es frekvenciákon). Ez egy alacsony-alacsony mőholdkövetésnek (low-low SST) felel meg. A két frekvencián meghatározott távolság nagyságrendekkel pontosabb a GPS mőholdak által meghatározott távolságnál (mm-es megbízhatóság) a megfelelı korrekciók figyelembe vételével. Célunk, hogy ezen pontosabb mérési eredményeket használjuk fel a számításainkhoz, de elsı megközelítésben (ezen tanulmány keretén belül) a bázisvonalra szorítkozunk.
1.3 Perturbációk A Földet elsı közelítésben tekinthetjük pontszerőnek, melynek tömege M, így egy gömbszimmetrikus erıteret képez (R = 6378137 méter). A potenciál értékének kiszámítása a következı képlettel lehetséges: GM g = 2 , melybıl látható, hogy ez a tömegközépponttól való távolság négyzetével r fordítottan arányos, ahol a G a gravitációs állandó. 2 Nm 2 14 Nm GM =3,986*10 ) (M=5,974*1024 kg, G=6,672*10-11 kg 2 kg A Föld gömbszimmetrikus erıterétıl való eltérés több tényezı függvényeként írható le: I. A Föld tömegeloszlása nem gömbszimmetrikus, melyet szférikus gömbfüggvénysorban fejezhetünk ki. l
GM ∞ R l ∑ ∑ (Clm ⋅ cos mλ + Slm ⋅ sin mλ )Plm (sinψ ) R l =0 r m=0 Ahol r, φ és λ geocentrikus polárkoordináták. Az l a szférikus harmonikus rendje, az m a fokszámát jelenti a függvénynek. U (r ,ϕ , λ ) =
4
II. III. IV. V.
VI.
VII.
A többi égitest perturbáló hatása, mint a Hold, a Nap és a többi bolygó (Luniszoláris perturbációk). A Föld lapultságának indirekt hatása, melybıl kifolyólag ez perturbációt okoz a Hold mozgásában, és emiatt eltolódik a Föld-Hold rendszer tömegközéppontja. Az általános relativitáselméletbıl adódó korrekció, amely a Földtıl, mint tömegponttól származó gyorsulással arányos, ahol az arányossági tényezı a Föld Schwarzschild-sugarának és a hold r távolságának hányadosa. Légköri fékezıdésbıl származó perturbációk. Ezen perturbációkat okozó erık két erı komponensre és egy nyomatékra bonthatók: a) Közegellenállási erı (fD), mely a haladási iránnyal ellentétes irányban hat. 1 f D = ρ ⋅ Vr2 ⋅ A ⋅ C D 2 Ahol: ρ a levegı sőrősége, A a holdnak a mozgás irányára merıleges keresztmetszete és CD a dimenziótlan közegellenállási együttható. b) Aerodinamikai felhajtóerı (fL) c) Forgatónyomaték (M) A sugárnyomásból származó perturbációk két részre oszthatók: 1) Direkt sugárnyomás: A Nap elektromágneses sugárzása a mesterséges holdra nyomóerıt fejt ki. A sugárzás részben elnyelıdik, részben visszatükrözıdik, és részben szóródik. A Föld árnyékában azonban nincs napsugárzás, azaz itt másként kell figyelembe venni ezt a hatást. Egy árnyékfüggvény bevezetése célravezetı, azonban nem teljes értékő figyelembe vétel. 2) Indirekt effektusok: a) A Földrıl visszaverıdött sugárzás sugárnyomása. Ez a Nap sugárzásának a Földrıl tükrözésszerően visszaverıdött, illetve diffúzan szórt része, valamint a Föld termikus sugárzása. b) A mesterséges hold anizotrop termikus emissziója. Azaz a hold hımérséklet-eloszlása nem egyenletes, aszimmetria lép fel a hold megvilágított és az éjszakai oldala között. c) A mesterséges hold rádiósugárzása. A kibocsátott sugárzás, mintegy hátralöki a holdat. Árapály perturbációk. A Hold és a Nap által a Földön okozott árapály olyan jelenségeket idéz elı, melyek nem hagyhatók figyelmen kívül a mesterséges holdak pályaszámításaiban. 1) Az árapály erık periodikus pulzációra késztetik a Földet (kinematikai effektus). 2) Az árapály erık a geopotenciál idıbeli változását okozzák, ami visszahat a mesterséges holdak mozgására is (dinamikai effektus). 3) Az árapály erık befolyásolják a Föld tengely körüli forgását, ez a pályszámításban használt koordináta rendszerekre van hatással (referencia rendszer effektus).
5
2 A közvetítıegyenlet levezetése Általános értelemben a munka, illetve az energia (T) a helynek (x) közvetett, a sebességnek ( x& ) közvetlen függvénye: T = T ( x, x& ) Mint tudjuk a munka az azt kialakító erırendszer (Q) hely szerinti integrálja: T = ∫ Qd x
Ebbıl az erıt (Q) deriválással nyerjük. Mivel T = T ( x, x& ) , így a teljes differenciál mind x, mind x& szerinti parciális tagokat tartalmaz. A munka teljes differenciálját a Lagrange-egyenlet adja meg: d ∂T ∂T =Q, − dt ∂ x& ∂ x
ahol a potenciál negatív elıjele fizikai konvenciókon alapszik (ellentétben a geodéziában használatos elıjellel, vö. 211).
2.1 Energia-integrál A kéttest-probléma feladata két tömegpont mozgásának vizsgálata, melynek során rájuk csak a Newton-féle kölcsönös gravitációs vonzóerı hat. A két test mozgása háromféle egyenlettel írható le: • impulzusmomentum-integrál • energia-integrál • Laplace-integrál Az energia-integrál segítségével egy egységnyi tömegő pont energia megmaradása a tömegpont kinetikus energiájával és a potenciális energiájának összegeként írható le: 1 &2 µ r − =h 2 r Ezen rendszer összenergiája (H) a tömegegységre vonatkoztatott energia (h) segítségével a következıképpen írható fel: H=
m1m2 h m1 + m2
Bontsuk fel az energiát (H) egy helyfüggı (V) és egy sebességfüggı (T) komponensre, amivel a Hamiltoni egyenlethez jutunk: (211)
H = T −V
A geodéziai gyakorlatban a sebességfüggı komponenshez (T) képest a helyfüggı komponenst (V) negatív elıjelőnek definiáljuk a klasszikus fizikai értelmezéssel szemben [vö. Lagrange egyenlet]. A sebesség függı komponenst kinetikus energiának, a helyfüggıt pedig potenciális energiának (potenciálnak) nevezzük.
6
A helyzeti energia (V) az erıtérben nem állandó, mivel a külsı (F) és belsı erık ( &x& ) hatására folyamatosan változik. ∂V = &x& − F ∂x
(212)
A kinetikus energia definíció szerint T =
1 2 x& , ennek sebesség szerinti deriváltja: 2
∂T = x& ∂ x&
Az összenergia H = H ( x, x& , t ) alakú, így a teljes differenciálját a következıképp kapjuk: dH ∂H d x ∂H d x& ∂H = + + ∂ x dt ∂ x& dt dt ∂t
A (211) egyenlet alapján: dH ∂V ∂T d x& ∂V =− x& + − dt ∂x ∂ x& dt ∂t
A fent bevezetett összefüggéseket felhasználva kapjuk, hogy ∂H d x& d x& ∂V =− x& + F x& + x& − , ∂t dt dt ∂t
ami a
d x& x& tagok kiesése után: dt ∂H ∂V = F x& − ∂t ∂t
Ezt idı szerint kiintegrálva kapjuk: H = ∫ F x&dt − ∫
∂V dt + E0 ∂t
(213)
Ahol E0 az integrálási konstans. Könnyen beláthatjuk, hogy: dx
∫ F x&dt = ∫ F dt dt = ∫ Fd x Ezzel a potenciálra rendezett a (211) egyenlet az alábbi lesz: V =T −H =
1 2 ∂V x& − ∫ F dx + ∫ dt − E0 2 ∂t
Ahol: 1 2 x& 2
a kinetikus energia
∫ F x&dt
nem konzervatív erık hatása (azaz külsı erık)
∂V
∫ ∂t dt
idıben változó tagok hatása
7
(214)
E0
energia állandó
A hely függvényében V = V (C nm , S nm , r , ϕ , λ ) tehát megállapíthatóak a C nm , S nm normalizált gömbfüggvény együtthatók. Az idıben változó tagok hatása a következı perturbáló erıktıl függ: V = Vrotating earth + Vlunar tide + Vsolar tide + Vplanetary tides + Vsolid earth tide + Vocean tide + Vatmospheric tide + Vocean loading + Vatmospheric loading + Vother mass redistirbutions Ahol: Vrotating earth
forgó Föld hatása
Vlunar tide
Hold árapálya
Vsolar tide
Nap árapálya
Vplanetary tides
bolygók hatása
Vsolid earth tide
merev Föld hatása
Vocean tide
óceánok hatása
Vatmospheric tide
atmoszféra hatása
Vocean loading
óceánok mozgása
Vatmospheric loading
atmoszféra mozgása
Vother mass redistributions
egyéb tömegvonzási tagok
Elsı közelítésben tekintsük csak a Föld forgásának hatását, és hanyagoljuk el a precessziót, nutációt és a pólusmozgást, hiszen ezen tagok nagysága 10-12 - 10-13-on nagyságrendő csak, míg a Föld forgási szögsebessége: ωe = 7,292115 × 10-5 rad/s
A földrajzi hosszúság tehát egyenletesen változik, λ = −ωt , ezért:
dλ = −ω dt A V = V (r ,ϕ , λ ) potenciálfüggvény [vö. gömbfüggvénysoros alak] teljes differenciálja
∂V ∂V dr ∂V dϕ ∂V dλ = + + ∂t ∂r dt ∂ϕ dt ∂λ dt így alakul: ∂V ∂V dλ ∂V ∂V = = −ω = −ω ∂t ∂λ dt ∂λ ∂α
(215)
A hely derékszögő vízszintes koordinátái: x1 = r ⋅ cos(ϕ ) ⋅ cos(α ) x2 = r ⋅ cos(ϕ ) ⋅ sin(α )
8
Ezek rektaszcenzió szerinti parciális deriváltjai ∂x1 = − r ⋅ cos(ϕ ) ⋅ sin(α ) = − x2 ∂α ∂x2 = r ⋅ cos(ϕ ) ⋅ cos(α ) = x1 ∂α Tehát a rektaszcenzió szerinti parciális derivált operátora:
dx dx ∂ ∂ dx1 ∂ dx2 = + = − x2 ⋅ 1 + x1 ⋅ 2 ∂α ∂x1 dα ∂x2 dα dα dα A potenciálfüggvény idı szerinti változása [lásd (215)] ∂V ∂V ∂V ∂V = −ω = −ω − x 2 + x1 ∂t ∂α ∂x1 ∂x 2
a fentiek alapján: ∂V ∂V ∂V = −ω x1 − x2 ∂t ∂ x ∂ x 2 1
(216)
A V potenciálfüggvényt bontsuk fel egy földi és egy egyéb potenciálra:
V = VFöld + δV A rendszerre ható összes erı (F): F = ∇V + F
Ahol: ∇V
a rendszerben ható belsı erık
F
a rendszerre ható külsı erık
A V = VFöld + δV képletbe behelyettesítve a (212) azonosságot kapjuk:
∂VFöld ∂δV + = &x& − F ∂x ∂x (216) alapján az idı szerinti változása a potenciálfüggvénynek:
∂V ∂δV = −ω x1 &x&2 − F2 − ∂t ∂x2 Mivel &x&i >> F i +
∂δV − x 2 &x&1 − F1 − ∂ x1
∂δV ez leegyszerősödik: ∂ xi
∂V ≈ ω ( x 2 &x&1 − x1 &x&2 ) ∂t
(217)
Amely a Föld egyenletes forgásából származó potenciál.
9
Kiintegrálva idı szerint:
dVω = ω (x 2 &x&1 − x1 &x&2 )dt Lássuk be, hogy
d (− x 1 x& 2 ) d ( x 2 x& 1 ) = − x& 1 x& 2 − x 1 &x&2 és = x& 1 x& 2 + x 2 &x&1 , ezért: dt dt
Vω = −ω (x1 x& 2 − x 2 x& 1 ) Összegzésül felírjuk egy test belsı energiájának helyfüggı komponensére származtatott közvetítıegyenletet: ∂V perturbáló 1 2 V = x& − ∫ F x&dt − ω ( x1 x& 2 − x 2 x& 1 ) + ∫ dt − E0 (218) 2 ∂t (Amely megfelel a CHAMP mőholddal végzett mérések közvetítıegyenletének.)
10
2.2 Low-Low Satellite to Satellite Tracking egyenlete 2.2.1 Low-low SST bázisvonalon alapuló megoldásának közvetítıegyenlete A kinetikus energiák különbsége a két mőhold között:
∆T = TB − T A =
1 1 2 x& B − x& A 2 2
2
melyet átrendezve kapjuk:
∆T =
1 [(x& B + x& A )(x& B − x& A )] 2
(221)
Ezen a ponton két folytatási lehetıség van, attól függıen, hogy melyik mőholdat tekintjük elsıdlegesnek: 1. x& B + x& A = x& B + x& B − x& B + x& A = 2 x& B − x& AB ekkor: 1.eset =
1 1 2 T x& AB (2 x& B − x& AB ) = − x& AB + x& B x& AB 2 2
2. x& B + x& A = x& B + x& A − x& A + x& A = 2 x& A + x& AB ekkor: 2.eset =
1 2 1 2 T x& AB (2 x& A + x& AB ) = x& AB + x& A x& AB 2 2
Ebbıl két a két mőholdra vonatkozó közvetítıegyenlet: V AB =
1 2 T A B B A A B B A x& AB + x& A x& AB − ∫ F B x& B − F A x& A dt − ω ( x12 x& 2 − x 2 x&12 − x1 x&12 + x12 x&1 ) − E0AB 2 t (222)
(
)
(Ez megfelel a GRACE bázisvonalon alapuló megoldás közvetítıegyenletének.) Megjegyzések:
∂V perturbáló
•
∫
•
A GRACE közvetítıegyenlete annyiban szolgáltathat pontosabb eredményt a
∂t
dt tag mindkét mőholdra nézve közel azonos, így elhanyagolható.
potenciálra (a gömbfüggvény együtthatókra), mint a CHAMP, amennyivel 1 1 2 2 T x& AB + x& A x& AB pontosabb x& -nél. 2 2
11
2.2.2 Low-low SST mikrohullámú távolságmérésen alapuló megoldásának közvetítıegyenlete A mőhold-mőholdkövetés gyakorlatilag két mőhold kölcsönös távolságmérése egymásra. A továbbiakban a mőholdak közötti pontos mikrohullámú távolságmérés alkalmazásának módját vezetjük le. Ehhez elıször értelmeznünk kell a távolságváltozás mérések fizikai tartalmát. Két mőhold közötti távolság a következıképpen írható fel:
ρ12 = x12 = e12 T x12 = x12 2 + y12 2 + z12 2 Az alsó ábra szerinti értelmezésben az egymástól konstans távolságra haladó mőholdak esetén (ui. v1 ≅ v 2 ) az e&12 és e12 egységvektorok merılegesek egymásra:
T
e12 e12 = 1
(lásd [7])
T e&12 e12 = 0 T Ezért a két mőhold közötti távolság változásából ( ρ&12 ) a e&12 x12 tag kiesik:
ρ&12 =
dρ12 T T T = e&12 x12 + e12 x& = e12 x& 12 dt
(223)
A potenciál (218) egyenlete egyszerősítı feltételezések mellett (atmoszferikus fékezıdés mentes pálya, statikus gravitációs tér) az alábbi alakra módosul:
V=
1 2 x& − E0 2
Ennek megváltozását fogjuk a továbbiakban koordinátarendszerében vizsgáljuk, melynek fıirányai: at: along track (pályamenti) ct: cross track (keresztirányú)
r: radial (sugárirányú)
elemezni.
Ezt
a
mőhold
X: a pozitív irány a keringés iránya Y: a pozitív irányt úgy választjuk, hogy jobbsodrású rendszert kapjunk ( e ct = e r × e at ) Z: a pozitív irány a Föld tömegközéppontja felé mutat
A potenciál pályamenti változása: (224)
d aV = x& T d a x&
12
A két egymáshoz közel keringı mőholdra (az ábrán 1 és 2 helyzet) beláthatjuk, hogy a pályamenti potenciál változás a pontbeli potenciálok különbségével közelíthetı, tehát d aV ≈ V2 − V1 Továbbá a pályamenti sebesség változás a mőholdak közötti sebességváltozás pálya irányú komponensével helyettesíthetı, ami pedig (223) alapján a távolságváltozással, tehát: T
d a x& ≈ e12 x& 12 = ρ&12
Ezek alapján a GRACE mőholdakra (A és B mőholdak) a (224) összefüggés az alábbi közelítı alakot nyeri: (225)
V AB = x& A ρ& AB
Megjegyezzük, hogy a levezetés során használt matematikai feltételezés jogosultsága, mely szerint a mőholdak közötti távolság kicsi, ρ AB → 0 , nem igazolható. A (225) egyenlet segítségével módunkban áll a távolságváltozás ( ρ& AB ) méréseinek felhasználására a közvetítı egyenletben. A továbbiakban a mőholdak közötti távolságmérés ( ρ AB ) értékeit is szeretnénk felhasználni, arra az esetre, hogyha ezekbıl a bias* értékét korrekten el tudjuk távolítani. Ezek alapján a két mőhold közötti távolság: x12 = ρ e12 ,
melynek teljes differenciálja: x&12 = ρ& e12 + ρ e&12
Ezt a differenciált négyzetre emelve a következı egyenlethez jutunk: 2 x&12 = ρ& 2 e12 e12 + 2 ρρ& e12 e&12 + ρ 2 e&12 e&12
(226)
További azonosságok: e12 × e12 = 0 e12 × e&12 = e ct
Ezek után lássuk be:
ρ −2 (( x12 × x&12 ) × e12 ) =
(227)
= ρ −2 ((ρ e12 × ρ& e12 + ρ e12 × ρ e&12 ) × e12 ) = ρ −2 ((ρ e12 × ρ& e12 + ρ e12 × ρ e&12 ) × e12 )
( (ρ
= ρ −2 ρ 2 e ct × e12 = ρ −2
2
)
e&12 )
= e&12 *
bias - mikrohullámú távolságmérés ismeretlen eltolás értéke (ami az észlelési eljárás következménye)
13
A (226) egyenletbe behelyettesítve a (227) egyenletet: 2 2 x& 12 = ρ& 2 + 2 ρρ&ρ −2 (( x12 × x&12 ) × e12 )e12 + ρ 2 ρ −2 ρ −2 (( x12 × x&12 ) × e12 )
(228)
Ahol:
( x 12
× x& 12
)=
ρ ρ& e ct
2 ρρ&ρ −2 (( x12 × x&12 )× e12 )e12 =
= 2 ρ&ρ −1 (( x12 × x&12 ) × e12 )e12 = 2 ρ&ρ −1 (ρρ& e ct × e12 )e12 = 2 ρ&ρ −1 (ρρ& e&12 )e12 = 2 ρ& 2 (e&12 e12 ) =0
A (228) azonosság a következı alakra egyszerősödik: 2 2 x& 12 = ρ& 2 + ρ −2 (( x12 × x&12 ) × e12 )
(229)
A (225) és a (229) összefüggéseket alkalmazva (222)-re a közvetítı egyenlet így alakul: V AB =
1 2 1 −2 (x AB × x& AB )2 + x& TA ρ& AB ρ& AB + ρ AB 2 2
(
)
A B B A A B B A − ∫ F B x& B − F A x& A dt − ω ( x 12 x& 2 − x 2 x& 12 − x 1 x& 12 + x 12 x& 1 ) − E 0AB
(2210)
t
(Ez megfelel a GRACE mikrohullámú távolság- és távolságváltozás mérésen alapuló megoldás közvetítıegyenletének.) Megjegyzések:
∂V perturbáló
•
∫
•
A GRACE közvetítıegyenlete annyiban szolgáltat pontosabb eredményt a
∂t
dt tag mindkét mőholdra nézve közel azonos, így elhanyagolható.
potenciálra (a gömbfüggvény együtthatókra), mint a CHAMP, amennyivel 1 2 1 2 1 −2 (x AB × x& AB )2 + x& TA ρ& AB pontosabb x& -nél. ρ& AB + ρ AB 2 2 2
14
3 Mérések feldolgozása 3.1 Adatok A GRACE mőholdak, valamint a mőholdakra vonatkozó mérési adatokat a Müncheni Mőszaki Egyetemen (TUM) a Csillagászati és Fizikai Geodéziai Intézetétıl (IAPG) kaptuk. A rendelkezésünkre bocsátott négyhónapnyi (2003. 07. 01. → 2003. 10. 31.), 5 másodpercenként rögzített mérési eredményekbıl 10 másodperces felbontásban meghatározott pályakoordináták napi felbontású fájljait ASCII formátumban tölthettük le.
3.2 Pályatípusok 3.2.1 Kinematikus pálya A nyers GPS mérési adatokból pályakoordináták számítása. Tulajdonságok: • Pusztán geometriai adatokból meghatározott pálya, mely független mérésnek tekinthetı a nehézségi erıtér meghatározása során. • Valóságos adatokat szolgáltat, így minden pontban mérési hibával terhelt.
3.2.2 Dinamikus pálya A GPS mérésekbıl csak egy pont koordinátáit, illetve sebességét használjuk fel kezdeti értékként, majd egy geopotenciális modell (Clm, Slm) felhasználásával kiintegráljuk a pályát. Tulajdonságok: • Pusztán a priori fizikai ismereteken nyugvó pálya • Kezdeti érték hibáival terhelt, melyek kiintegrálódnak, azaz halmozódnak. • Modellhibákkal is terhelt.
15
3.2.3 Féldinamikus pálya (reduce dynamic orbit) A kinematikus és a dinamikus pálya ötvözése, amikor a kezdeti érték feltétel mellett kényszerfeltételeket is alkalmazunk. Bizonyos feltételek mellett, pedig mindig új kezdeti érték feltételt határoz meg. Tulajdonságok: • Szakaszosan dinamikus, tehát sima pálya. • Ötvözi a dinamikus és a kinematikus pályák elınyeit; közel valós pálya és sima (mérési hibáktól mentes).
3.3 Programok A végsı eredmények elérése érdekében a programozás lépéseit részekre bontottuk az áttekinthetıség és a könnyebb kezelés érdekében. 1) Sebességszámítás 2) Adatbeolvasás és feldolgozás 3) Energia-integrál kiszámítása 4) GRACE-EGM potenciálkülönbség meghatározása a drift levonásához 5) A hosszú hullámhosszú tag meghatározása és levonása 6) Kiegyenlítés és a kapott eredmények elemzése
16
3.3.1 Sebességszámítás A kinematikus pályák sajátossága, hogy pusztán pozíció adatokat szolgáltatnak, ezért a kinetikus energiában szereplı sebesség más úton számolandó. Ezt 7-rendő NewtonGregory interpolációs eljárással analitikusan határoztuk meg egy mára gyakorlatban bevált program használatával.
3.3.2 Adatbeolvasás és feldolgozás A program feladata a dinamikus pozíció és sebesség, valamint a kinematikus pozíció adatok beolvasása, illetve ezek feldolgozása. Ennek a kezdeti feltételei azon két idıpont, melyek között szeretnénk elvégezni a fent említett mőveleteket (2003. 07. 01. → 2003. 10. 31.). A program alapján a kezdeti évet és az évben a nap számát (1-365), illetıleg az utolsó évet és az évben a nap számát (1-365), valamint a lépésközt kellett ehhez megadnunk. A KBR adatok beolvasását követıen a kinematikus pozíció adatokat szőrésnek vetettük alá, mert a KBR adatok felbontása, 5 s, eltér a pályaadatok felbontásától, 10 s (down sampling). Mivel nem minden adatát használtuk fel a beolvasott állományoknak az általunk elvégzendı számításokhoz, így egy-egy fájl beolvasása után ezeket is töröltük, segítve ezzel a számítógép memóriájának jobb kihasználását. Az adatok beolvasását követı mőveletek: • A tényleges adatfeldolgozás elsı lépéseként a GPS idırendszerébıl (gpstime) létrehoztunk az általunk vizsgált idıintervallumra egy folytonos idısort (tK), mely 5 másodperces felbontású, és 0 másodperccel indul (2003. 07. 01. 0: 00: 00). • A két mőhold között mért távolságok és távolságváltozások javítása a futási idı korrekcióval és az antenna fáziscentrum külpontosság korrekciójával. • Azon helyek megkeresése az adatsorban, ahol az adatok minıségét jelzı index fázisugrást mutatott, illetıleg azon helyeké, ahol az adatszolgáltatók rossz adatminıséget jeleztek (tehát quality flag = 1). Ezeket a mérési eredményeket töröltük a memóriából. • A KBR adatok folyamatos, 5 másodperces felbontású tK idıskálához való igazítása. • Újabb adatellenırzés, amely a KBR adatokban található adathiányos helyek indexeinek felkeresését végzi. (Figyelemmel kellett arra lenni, hogy az indexsor az egész adatsort lefedje!) • A KBR adatok szakaszokra bontása a hirtelen bias ugrás helyeken, mely a fázisugrások helyeit és a KBR adatokban lévı üres helyeket veszi figyelembe. (Több száz kilométeres bias is volt!)
•
A 10 másodperces felbontású geometriai pálya beolvasása, mely a GPS mőholdak által mért adatok feldolgozott pályaadatait tartalmazza a GRACE A és GRACE B mőholdakra.
17
• •
•
•
•
•
•
•
• • •
A dinamikus pályaadatokból szintén definiálunk egy folytonos idıskálát (tA), mely 10 másodperces felbontású. Megkeressük ebben az idısorban azokat a helyeket, ahol nincs adat (NaN), mely esetünkben egy összefüggı 8640 db-ból álló adatsor volt, ami azt jelenti, hogy egy teljes nap hiányzik a dinamikus pályaadatokból. A dinamikus idısort (tA) ezután két részre bontottuk ennél az adatszakadásnál (tA1 és tA2), valamint a KBR adatsorból meghatározott kinematikus idısort (tK) is ezen határok mentén bontottuk két részre. (tA 10 s-os, tK 5 s-os felbontású!) Egy vizsgálatot végzünk ezután, melynek a feladata, hogy ellenırizze a két GRACE mőhold idıvektorának egyezıségét. Ez a feltétel a mi esetünkben teljesült. Megkeressük azokat a koordinátákat a GRACE A és GRACE B mőholdakra vonatkozólag a dinamikus pályaadatokból, amelyekre számszerő értéket kaptunk. A két részre osztott adatsornál alkalmazunk egy interpolációs eljárást, mely a kinematikus pályaadatokat felhasználva ugyanazon idıpontokra (tA1-rıl tK1-re valamint tA2-rıl tK2-re) interpolál spline függvény segítségével dinamikus távolságot (r0) a mőholdak között. (Egy adott idıpont csak egyszer fordulhatott elı, és ehhez csak egy távolság tartozhat!) A két mőhold közötti interpolált és eredeti féldinamikus távolságok összehasonlítása. Esetünkben a két függvény karakterisztikáját nézve jól megegyeztek, tehát az interpoláció hibái elhanyagolhatóak.
Adatszőrés, mely az adathiányos helyek felkeresését végzi a dinamikus pályaadatokban (lásd fenti ábra). Ezeken a helyeken tördeljük az adatsorokat a bias becsléséhez. A biasok becslése az egyes szakaszokra számtani közepet használva értékének meghatározására. (bıvebben lásd 3.2.2.1 pont) A bias értékek megjelenítése (lásd 3.3.2.2 A bias értékek nagysága) A szükséges adatok elmentése (távolságok, bias, pozíciók, sebességek, indexek).
18
3.3.2.1 A biasok becslése A bias meghatározás lényege, hogy a GRACE A és GRACE B mőholdak között dinamikus pályaadatokból számítható távolság, illetve a KBR mérések eredményei hasonló szinusz-görbéket írnak le, azon különbséggel, hogy a dinamikus pályából számított távolság folyamatosan megbízható eredményeket szolgáltat, viszont csak ±5 cm pontosan. A KBR mérésekben az egész rész meghatározása nagy bizonytalanságokat rejt magában, mivel elıfordul negatív távolság is, viszont a pontossága egy-egy szakaszon mm-es nagyságrendő, tehát nagyságrendekkel pontosabb, mint a dinamikus pályaadatokból számítható. Tehát van két közel szinuszos jellegő függvényünk; az egyik számított (GPS high-low SST mérésekbıl geopotenciál modell segítségével számított), a másik közvetlen mérési eredményekbıl származik. (Megjegyzés: A dinamikus pályát is mérési eredményekbıl számítottuk, amelyet az eddigi ismereteink alapján, geopotenciál modell alapján megerısítünk. A feladatunk, hogy meghatározzuk mennyi az eltolás a két függvény között egyegy szakaszon, és a mérési eredményt idetolva, ezek után ezt használjuk fel a további számításokhoz, mely tartalmazza a mérési zajokat, azaz a számunkra értékes információt. A feladat megoldását nehezíti, hogy az eltolás értéke nem állandó az idıtartam során, hanem ismeretlen helyeken (idıpontokban), ismeretlen mértékben változhat. Elsı lépésként meghatároztuk azokat a helyeket, ahol az eltolás értéke látványosan – akár vizuálisan is – észlelhetı. A bias meghatározási feladatot ezeken a helyeken értelemszerően külön választottuk, ezzel a folytonos távolságméréseket egymástól független szakaszokká bontva. A négy hónap adatai összesen 188 darab független szakaszt jelentettek. A probléma a kis mértékő bias változásoknál volt, ugyanis ezeket a helyeket nem lehetett észrevenni, és a bias értékében bekövetkezett ugrás hibaként bent marad a további számításokban. A legjobb módszer megállapítása érdekében többféle interpolációs eljárással próbálkoztunk. A 8 pontra illesztett 7-ed fokú polinom volt az elsı kísérletünk, mely a két függvény azonos idıbeli pontjaira egy polinomot illeszt, és ebbıl határoz meg egy közbensı pontra egy eltolás értéket. Ez a módszer a pontbeli bias (eltolás) meghatározásához csak a környezı mérések értékeit veszi figyelembe, függetlenül a teljes idısor hosszától és karakterisztikájától. Ezzel az adatból ki nem szőrhetı apró bias ugrások hibahatását szerettük volna kisebb, lokális területre korlátozni. Esetünkben ez azonban nem vált be, a lokális meghatározási eljárás miatt ugyanis a becsült bias értékek az eredeti függvénnyel megegyezı fázisban 150 méteres szinusz hullámú zajt generáltak. Ezt a megoldást elvetettük. A következıkben magasabb fokszámú polinommal próbáltuk meghatározni az eltolódás mértékét, melytıl a szinuszos jellegő függvény amplitúdójának csökkenését vártuk. 16 pontra illesztett 15-öd fokú polinommal a tesztterületre 300 méter nagyságú szinuszos hibát kaptunk.
19
A feltételezésünket megfordítva ezek után 3 pontra illesztett 2-od fokú polinommal próbálkoztunk, melynek eredményeként 50 méteres szinuszos hullám maradt meg mérési zajként, végül 2 pontra illesztett polinom (egyenes) esetén is 5 centiméter amplitúdójú szinuszos jellegő függvényt kaptunk, melynek alkalmazását nem fogadhattuk el. A nagyon pontos KBR mérési eredményeket csak abban az esetben tudjuk felhasználni, ha a bias meghatározása is hasonló a mérési eredmények pontosságához, azaz mm pontosságot kellett elérnünk. Ezek után – látva azt, hogy a függvények között az alacsonyabb fokú függvények adnak jobb megoldást – döntöttünk az egyszerő számtani közép képzése mellett, mely az adott szakaszokra vonatkozik, és nem csak a függvény néhány kiválasztott pontjára. Ezzel a megoldással a becsült biasról biztosan tudjuk, hogy a kisebb ugrásait nem veszi észre, értékei torzítottak. Pontosságának mértéke így megbecsülhetetlen, viszont mesterségesen szinuszos hibát nem vittünk bele. A megoldás pontosságáról visszacsatolás nincs, konklúziót csak a gömbharmonikus együtthatók szintjén adhatunk. Azon a szinten viszont már egyéb hibahatásokkal együtt érzékelhetjük csak a bias meghatározás hibáját, ezért egy biastól független eljárás használatát is szükségesnek tartjuk. Végeredményként tehát megkaptuk a KBR mérések nyújtotta nagyobb pontosságú távolságokat melyekben már nincsenek ugrások, és a továbbiakban ezen mérési eredményeket használhatjuk fel a számításokhoz.
3.3.2.2 A bias értékek nagysága A bal oldali ábrán kék vonallal látható a KBR mérések eredményeinek feldolgozásával kapott biasmentes távolság. Fekete vonallal a KBR mérések eredményei, és piros vonallal a bias. A jobb oldali ábrán a már eltolt szignál látható kinagyítva, mely karakterisztikáját tekintve pedig szinuszos jellegő. Jól megfigyelhetı, hogy a két mőhold közötti távolság 220 km körül változik,
3.3.3 Energia-integrál kiszámítása A további számításokhoz szükséges volt meghatároznunk a potenciálzavart, a tömegvonzási potenciált, a normál potenciált, a centrifugális potenciált és a kinetikus energiát. A program lépései: • Dinamikus pályaadatok beolvasása (10 másodperces felbontású). • Mivel a gyorsulásmérı adatait nem dolgoztuk még fel, ezért ezeket az adatokat 0-val tesszük egyenlıvé. • Az adatbeolvasás és feldolgozás során elmentett fájl beolvasása.
20
•
• • •
•
A KBR adatok igazítása a 10 másodperces idıfelosztáshoz (down sampling). Valamint a megfelelı távolságok és sebességek kiszámítása súlyozással. ai −1 + 2ai + ai +1 4 A szakaszok elején, illetve végén természetesen nem lehetséges ez a megoldás, így ott a becslést az arra az idıpontra vonatkozó adat kétszeres súlyával, és a szakasz elején az ezt követı 5 másodperces adat egyszeres súlyával, illetve a szakasz végén az ezt megelızı 5 másodperces adat egyszeres súlyával számítjuk. 2ai + ai +1 ai −1 + 2ai illetve 3 3 A kinematikus sebesség, pozíció és egy ezek megbízhatóságát minısítı bináris vektor (quality flag) számítása. A rossz minıségő adatok helyének (vö. quality flag), illetve az adatsorokban található NaN értékek indexeinek a felkeresése. Ábra nyomtatása annak érdekében, hogy ellenırizzük a GRACE pálya adatok és a KBR mérések idıpontjainak egyezıségét. Egy teljes napi KBR adat hiányzott, ezért itt meg kellett szakítani az adott szakaszt (208. nap hiányzott → 2003. július 27.). Energia-integrál segítségével szakaszokra bontva a potenciálzavar, a tömegvonzási potenciál, a normál potenciál, a centrifugális energia és a kinetikus energia kiszámítása és mentése.
3.3.3.1 A munka jellegő mennyiségek nagyságrendjei Jelen esetben a két mőhold közötti energiakülönbségek nagyságát tőntettük fel minden ábrán!
• • • •
Kinetikus energia Centrifugális potenciál Normál potenciál Potenciálzavar
σ = ±8100 m2/s2 σ = ±2777 m2/s2 σ = ±16,7 m2/s2 σ = ±17,2 m2/s2
21
3.3.3.2 Nem konzervatív erık nagysága Jelen esetben a két mőhold közötti energiakülönbségek nagyságát tőntettük fel minden ábrán!
• • • • •
A Nap és a Hold által okozott direkt árapály becslése Merev Föld árapálya Óceán árapálya Pólusmozgás Egyéb égitestek által okozott direkt árapály
22
σ = ±0,1 m2/s2 σ = ±0,02 m2/s2 σ = ±0,01 m2/s2 σ = ±0,001 m2/s2 σ = ±0,1 m2/s2
3.3.4 GRACE-EGM potenciálkülönbség meghatározása a drift levonásához A program feladata a GRACE mőholdak mérési eredményeibıl, és az EGM96-os modellbıl számolt pályamenti potenciálértékek különbségének meghatározása a drift levonásához. A program lépései: • Az EGM96 geopotenciál modell betöltése. • Egy rácsháló definiálása a gömbön, amely értékeket 0,25°-ként vettünk fel φr szélesség és λr hosszúság irányokban egyaránt.
•
•
Az elızıekben kiszámított adatok beolvasása után a derékszögő koordinátákat átszámítottuk gömbi koordinátákká, figyelembe véve azt, hogy nem radiánban, hanem tizedfokban kívánjuk felhasználni. (Ezt minden egyes szakaszra külön elvégeztük.) Ezután az egész idıtartamra (4 hónap) beolvastuk a GRACE A mőhold geocentrikus távolságait egyetlen vektorba, középértéket képeztünk ebbıl (rA), majd kiszámítottuk a GRACE A mőhold átlagos magasságát az általunk meghatározott 6378137 méter sugarú gömb felett (hm). Ezzel együtt meghatároztuk a szakaszokra különkülön, minden egyes felvett idıpontra az
23
• •
•
•
átlagtól való eltérését a geocentrikus távolságnak (dhA). Ugyanezen mőveleteket elvégeztük a GRACE B mőholdak mérési eredményire is. Végül a mindkét mőhold adatait felhasználva számítottunk egy közepes geocentrikus távolságot ( r ), és határoztuk meg az eltéréseket minden egyes távolságra. A következı lépésben gömbharmonikus szintézissel meghatároztuk az EGM96 modellbıl a potenciálértékeket a raszter φr, λr, r koordinátájú helyeire a gömbfüggvény együtthatók segítségével (2-360 fokig). Ezután az r sugarú gömb felszínén interpolációval határoztunk meg az adott φ, λ helyre potenciál értéket (még a gömb felszínén), végül magassági értelemben a Taylor-sor segítségével a gömb felszínérıl a megfelelı magasságba korrigáltunk. 2 ∂V (r ,ϕ , λ ) (r − r ) + 1 ∂ V (r ,2ϕ , λ ) (r − r )2 V (r ,ϕ , λ ) = V (r ,ϕ , λ ) + ∂r 2 ∂r Ezt a mőveletet minden egyes szakaszra, mind a két mőholdra külön-külön elvégeztük, és mentettük ezen számított potenciálértékeket. Végezetül ábrán szemlélhetjük meg a számított potenciálkülönbség értékeket. Minden esetben a vízszintes tengelyen a földrajzi hosszúság szerepel 1441 pontra felosztva, a földrajzi hosszúság pedig a függıleges tengelyen, mely 721 részre van felbontva.
24
3.3.5 A hosszú hullámhosszú tag (drift) meghatározása és levonása Ebben a lépésben számítottuk ki a drift értékét, majd vontuk le a potenciálzavar értékébıl. Végezetül egy idısorba illesztettük az összes változót.
A számítás lépései: • A fájlok neveinek definiálását itt is meg kellett oldani, hogy a MATLAB ezek után be tudja olvasni a szakaszonként elmentett adatsorokat. • Itt is szükséges volt, hogy az elmentett derékszögő koordinátákat átváltsuk gömbi koordinátákká, így ez lett a következı lépés. • Az adatok szelektálása. Amennyiben a driftben törés lenne, akkor evvel a lépéssel egyrészt két részre bonthatjuk az adatsort, és külön vizsgáljuk a driftjüket.
25
•
Esetünkben a hosszú hullámú tag mindösszesen néhány m2/s2 körül váltakozott, amelyben alig észrevehetı törések maradtak csak (ezeket nem lehetett tovább szőrni, mert olyan kicsik voltak – néhány tized m2/s2). A GRACE A illetve GRACE B mőholdakra meghatározott potenciál értékekbıl számítottunk potenciálkülönbséget. Valamint a rájuk ható perturbáló erık különbségét is képeztük, melynek nagysága szintén elhanyagolható lett volna (néhány tized m2/s2), hiszen a két mőholdra közel azonosak ezek a hatások.
•
Ez az ábra jól szemlélteti azt, hogy a két mőhold potenciálkülönbsége mellett (néhány 10 m2/s2) a perturbáló erık, melyekkel a bemeneti adatokat javítottuk meg, szinte elhanyagolható mértékőek. Végül a különbözı változókat egy-egy vektorba győjtjük és egy fájlba elmentjük ezeket. A származtatott potenciálkülönbség értéket a két mőhold pályájának a magasságában, pont a fél úton értelmezünk [lásd következı ábra].
3.3.6 A kiegyenlítés és a kapott eredmények elemzése Ezek alapján mindkét mőholdra adott idıpontokra geocentrikus helyvektorokat határoztunk meg illetıleg a mőholdak közötti vektort (s), melyek több tényezı függvényei:
• • •
simuló pályaelemek, nehézségi erıtér véges számú gömbfüggvény együtthatója, álláspont helyvektora. 26
s = s(Ω0, ω0, i0, a0, e0, T0; Jn, Clm, Slm; t, r)
A hosszabb idın keresztül végzett távolságmérés és pozíció meghatározások eredményeként felállítható javítási egyenlet kiegyenlítését többféle úton végezhetjük el. Mi a hagyományos legkisebb négyzetek módszere szerinti kiegyenlítés mellett döntöttünk, melynek során 35 napos idıtartam (2003.07.01.-2005.08.05., a 2003.07.28. nap kihagyása mellett) adataiból 12 fokig és rendig határoztuk meg a gömbfüggvény együtthatókat. (Ez az adat mennyiség és fokszám felelt meg egy PC kapacitásának.) A potenciálkülönbség a következıképpen írható fel:
∆V AB
GM = R
R ∑ l = 2 rB ∞
GM − R
l +1
l
∑ (C
R ∑ l =2 rA ∞
lm
⋅ cos mλ B + S lm ⋅ sin mλ B )Plm (sinψ B )
m =0 l +1
l
∑ (C
lm
⋅ cos mλ A + Slm ⋅ sin mλ A )Plm (sinψ A )
m =0
A legkisebb négyzetek módszere szerinti kiegyenlítés közvetítı egyenletekkel történı megoldását (II. kiegyenlítési csoport) felhasználva határoztuk meg a Clm, Slm gömbfüggvény együtthatókat. l + v = A⋅ x
ahol a javításokat minimalizáltuk:
v T v = min
Esetünkben az egyenlet tagjai: l +1 ∂∆V AB GM R cos(mλ ) P (sinψ ) − GM B lm B R rB R A = ∂C lm = l +1 GM ∂∆V AB GM R ∂ S lm R r sin(mλ B ) Plm (sinψ B ) − R B l = [∆V AB ]
l +1 R cos(mλ A ) Plm (sinψ A ) rA l +1 R sin( mλ A ) Plm (sinψ A ) rA
C lm x= S lm
A kiegyenlítéssel kapott együtthatókat a továbbiakban az alábbi elrendezés szerint mutatjuk. C00 C 10 C 20 C30 . . C n0
S11
S 21
S 31
C11
S 22
S 32
C 21 C 22
S 33
C31
C32
C33
. .
. .
. .
.
.
.
. . S n1 . . . . . . . . . . . . . . S nn . . C nn
A mellékelt ábrákon az együtthatók logaritmikus beosztású skálán láthatók. Az ábrák jobb oldalán a színskála mutatja az együtthatók nagyságrendjét.
27
Mint látható, a GRACE mőholdak méréseibıl kapott gömbfüggvény együtthatók (lásd bal oldali ábra) jó egyezést mutatnak az EGM96-os modell együtthatóival (lásd jobb oldali ábra).
Az együtthatók eltéréseit jobban szemlélhetjük a köztük lévı differencia mátrix segítségével, amelyet szintén logaritmikus beosztású ábrán, színekkel jelöltünk. Felhívnánk a figyelmet a lenti és a fenti ábrák eltérı színskála használatára.
A gömbfüggvény együtthatók felhasználásával elıállított geoid modellünk (lásd bal oldali ábra) is nagy hasonlóságot mutat az EGM96-os modellel (lásd jobb oldali ábra).
28
3.4 Eredmények A kinematikus pályaadatok feldolgozásával meghatároztuk a két mőholdra ható konzervatív és nem konzervatív erık keltette potenciálkülönbségeket tagonként. Majd a feldolgozás végén a drift levonása után, a kiegyenlítéssel 12 fokig és rendig meghatároztuk a gömbfüggvény együtthatókat, tehát elıállítottunk egy geopotenciális modellt (geoidmodellt). A kapott eredmények pontossága még elmarad az eddigi modellek pontosságától, mivel nem használtuk fel a mikrohullámú távolságméréseket a két mőhold között, amelyekkel ez biztosítható lenne (ez várhatóan a diplomamunka keretén belül kerül feldolgozásra). Azonban ez az eredmény összhangban van a nemzetközi tapasztalattal: a GRACE mőholdak pályáit nem high-low SST feldolgozásra, vagy bázisvonalon alapuló megoldásra találták ki. A feladatra sokkal hatékonyabbnak bizonyult a CHAMP pálya. A teljes, KBR mérési adatokat is felhasználó megoldáshoz ez az elsı, azzal matematikailag szinte megegyezı megoldása (vö 225 és 2210 összefüggések) egy elkerülhetetlen elsı lépés az adatok minıségének megismeréséhez. Modellünk tehát így is nagy hasonlóságot mutat az eddigi modellekkel, és elérte az elızetes elvárásoknak megfelelı pontosságot.
29
4 Mellékletek 4.1 Rövidítések CHAMP ESA GFZ GOCE GRACE IAPG KBR MATLAB NASA SLR SST TUM
CHAllenging Mini-satellite Payload European Space Agency GeoForschungZentrum Gravity field and steady-state Ocean Circulation Explorer Gravity Recovery And Climate Experiment Institut für Astronomische und Physikalische Geodäsie K Band Ranging MATrix LABoratory National Aeronautics and Space Administration Satellite Laser Ranging Satellite to Satellite Tracking Technische University München
4.2 Irodalomjegyzék [1]
Dr. Biró Péter: Felsıgeodézia (BSc) (Elektronikus jegyzet) Budapest, 2004. 126 p. [2] Case, Kelley – Kruizinga, Gerhard L. H. – Wu, Sien-Chong: GRACE Level 1B Data Product User Handbook Jet Propulsion Laboratory, California Institute of Technology JPL D-22027, 2004. Január 21. [3] Detrekıi Ákos: Kiegyenlítı számítások Tankönyvkiadó Vállalat, 1991. [4] Érdi Bálint: Mesterséges holdak mozgása Nemzeti Tankönyvkiadó, 1993. [5] Földváry Lóránt: A 2000-es évek elsı évtizede: a gravimetriai mőholdak korszaka Magyar Geofizika 45.évf. 4. sz. 1-7. o. 2005. [6] Földváry Lóránt – Fukuda, Yoichi: On the effects of the atmospheric correction of the GRACE measurements for studies of oceanography Periodica Polytechnika ser. Civ. Eng. Vol. 46, No. 2, pp 185-198. 2002. [7] Jekeli, Christopher: The Determination of gravitational potential differences from Satellite-to-Satellite Tracking Celestial Mechanics and Dynamical Astronomy 75: 85-101, 1999. [8] Mayr, Theresa: Schwerefeldanalyse der Satellitenmission GRACE unter Verwendung des Energieintegrals Diplomarbeit, TUM, 2005. [9] Proposal for a German Priority Research Program: Mass Transport and Mass Distribution in the Earth System Contribution of the New Generation of Satellite Gravity and Altimetry Missions to Geosciences GOCE-Projectbüro Deutschland, Technische Universität München GeoForschungZentrum Postdam, 2005. Január [10] Sneeuw, Nico: A semi-analytical approach to gravity field analysis from satellite observations Deutsche Geodätische Kommission bei der Bayerischen Akademie der Wissenschaften, Dissertationen, Heft Nr. 527 München, 2000. 30
4.3 Feldolgozott adatok formátuma KZ21_03182.KIN
(pályaadatok fáljai)
GRACE zero-difference kinematic orbit, day 03182 06-JUL-04 14:53 -----------------------------------------------------------------------------------------------LOCAL GEODETIC DATUM: IGS00 EPOCH: 2003-07-01 00:00:00 RMS= 0.003259 STATION NAME WEEK SECONDS X (M) Y (M) Z (M)
F
GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD GRACE-A-POD
K K K K K K K K K K K K K K K K K K K K K K K K K
1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225 1225
172800. 172810. 172820. 172830. 172840. 172850. 172860. 172870. 172880. 172890. 172900. 172910. 172920. 172930. 172940. 172950. 172960. 172970. 172980. 172990. 173000. 173010. 173020. 173030. 173040.
3972911.4254 3914012.0194 3854609.2717 3794710.3663 3734322.5876 3673453.2788 3612109.9265 3550299.9864 3488031.0365 3425310.6703 3362146.5528 3298546.4282 3234518.0633 3170069.3175 3105208.0535 3039942.2286 2974279.8338 2908228.9282 2841797.6029 2774994.0210 2707826.3406 2640302.8141 2572431.7610 2504221.4917 2435680.3804
KBR1B_2003-07-01_X_00.asc PRODUCER AGENCY : PRODUCER INSTITUTION : FILE TYPE ipKBR1BF : FILE FORMAT 0=BINARY 1=ASCII : NUMBER OF HEADER RECORDS : SOFTWARE VERSION : SOFTWARE LINK TIME : REFERENCE DOCUMENTATION : SATELLITE NAME : SENSOR NAME : TIME EPOCH (GPS TIME) : TIME FIRST OBS(SEC PAST EPOCH): TIME LAST OBS(SEC PAST EPOCH) : NUMBER OF DATA RECORDS : PRODUCT CREATE START TIME(UTC):
1067085.7285 1050469.9705 1033811.8020 1017114.0543 1000379.5169 983610.9500 966811.0988 949982.7184 933128.5554 916251.3565 899353.8145 882438.6581 865508.5998 848566.3369 831614.5435 814655.8960 797693.0650 780728.6959 763765.4448 746805.9346 729852.7690 712908.5500 695975.8901 679057.3442 662155.4810
5492653.7122 5537938.0825 5582539.6918 5626453.0569 5669672.7344 5712193.4643 5754010.0963 5795117.4988 5835510.6356 5875184.5423 5914134.4392 5952355.5440 5989843.1509 6026592.7980 6062599.9055 6097860.1169 6132369.0971 6166122.6338 6199116.6791 6231347.1914 6262810.1797 6293501.8061 6323418.4717 6352556.3601 6380912.0046
(KBR mérések fáljai)
NASA JPL 7 1 47 @(#) KBR_compress.c 1.69 05/29/03 @(#) 2003-09-11 14:15:51 glk j2 GRACE Level 1 Software Handbook GRACE A+B IPU 1+1 2000-01-01 12:00:00 110289600.000000 (2003-07-01 00:00:00.00) 110375995.000000 (2003-07-01 23:59:55.00) 17141 2003-09-21 14:43:50 by l0tol1
31
PRODUCT CREATE END TIME(UTC) : 2003-09-21 14:44:00 by l0tol1 FILESIZE (BYTES) : 4170092 FILENAME : KBR1B_2003-07-01_X_00.asc PROCESS LEVEL (1A OR 1B) : 1B INPUT FILE NAME : KBR1A_A_0<-KBR1A_2003-07-01_A_00.dat INPUT FILE TIME TAG (UTC) : KBR1A_A_0<-2003-09-04 06:55:15 by l0tol1 INPUT FILE NAME : KBR1A_B_0<-KBR1A_2003-07-01_B_00.dat INPUT FILE TIME TAG (UTC) : KBR1A_B_0<-2003-09-20 19:57:37 by l0tol1 INPUT FILE NAME : KBR1A_A_1<-KBR1B_2003-07-01_A_00.dat INPUT FILE TIME TAG (UTC) : KBR1A_A_1<-2003-09-21 14:43:30 by l0tol1 INPUT FILE SOFTWARE VERSION : KBR1A_A_1<-@(#) KBR_debreak.c 1.46 05/29/0 INPUT FILE LINKTIME TAG : KBR1A_A_1<-@(#) 2003-09-11 14:15:51 glk j2 INPUT FILE NAME : KBR1A_B_1<-KBR1B_2003-07-01_B_00.dat INPUT FILE TIME TAG (UTC) : KBR1A_B_1<-2003-09-21 14:43:37 by l0tol1 INPUT FILE SOFTWARE VERSION : KBR1A_B_1<-@(#) KBR_debreak.c 1.46 05/29/0 INPUT FILE LINKTIME TAG : KBR1A_B_1<-@(#) 2003-09-11 14:15:51 glk j2 INPUT FILE NAME : KBR1A_A_2<-KBR1B_2003-07-01_A_00.dat.ord INPUT FILE TIME TAG (UTC) : KBR1A_A_2<-2003-09-21 14:43:43 by l0tol1 INPUT FILE SOFTWARE VERSION : KBR1A_A_2<-@(#) KBR_order.c 1.14 05/29/03 INPUT FILE LINKTIME TAG : KBR1A_A_2<-@(#) 2003-09-11 14:15:52 glk j2 INPUT FILE NAME : KBR1A_B_2<-KBR1B_2003-07-01_A_00.dat.ord INPUT FILE TIME TAG (UTC) : KBR1A_B_2<-2003-09-21 14:43:50 by l0tol1 INPUT FILE SOFTWARE VERSION : KBR1A_B_2<-@(#) KBR_order.c 1.14 05/29/03 INPUT FILE LINKTIME TAG : KBR1A_B_2<-@(#) 2003-09-11 14:15:52 glk j2 INPUT FILE NAME : PCI1A_A<-PCI1A_2003-07-01_A_00.dat INPUT FILE TIME TAG (UTC) : PCI1A_A<-2003-09-21 14:33:34 by l0tol1 INPUT FILE NAME : PCI1A_B<-PCI1A_2003-07-01_B_00.dat INPUT FILE TIME TAG (UTC) : PCI1A_B<-2003-09-21 14:42:00 by l0tol1 INPUT FILE NAME : USO1B_A<-USO1B_2003-07-01_A_00.dat INPUT FILE TIME TAG (UTC) : USO1B_A<-2003-09-21 14:35:47 by l0tol1 INPUT FILE NAME : USO1B_B<-USO1B_2003-07-01_B_00.dat INPUT FILE TIME TAG (UTC) : USO1B_B<-2003-09-21 14:40:34 by l0tol1 END OF HEADER 110289600 -111586.0534618485 0.2742853150043325 0.002283814079108859 282354.504924697 -9.806928322154786e-05 8.385911015226631e-07 5.448293973907194e-09 2.945149543059169 7.372721961501486e-10 5.995844449961493e-08 679 688 698 703 00000000 110289605 -111584.6534700724 0.2857151143034363 0.002288216327423613 282354.5049294686 -9.387131262383961e-05 8.403307899273224e-07 2.901133028471517e-10 2.945149261758897 -5.193350163129979e-08 4.956176960083296e-09 671 674 698 703 00000000 110289610 -111583.1962753954 0.2971657718430906 0.002291899989830055 282354.504949363 -8.966607202534021e-05 8.417552599262692e-07 2.796790400725994e-10 2.945148656994503 -1.560402333053762e-07 3.606625498202343e-08 671 674 683 685 00000000 110289615 -111581.6817819598 0.3086348157738243 0.002295693255697416 282354.5049698481 -8.545386866933774e-05 8.43114810276783e-07 2.634615019637271e-10 2.945148631125216 8.573512371182807e-08 2.178580141919452e-08 682 688 683 685 00000000 110289620 -111580.1098987928 0.3201208834721528 0.002298661188666356 282354.5049789351 -8.123504287329455e-05 8.444032354555723e-07 2.498891976047794e-10 2.945148528572038 -7.820304095620402e-08 2.040771194395257e-08 682 688 697 706 00000000
32
4.4 A mőholdak pillanatnyi helyzete A http://www.csr.utexas.edu/grace/ honlapján 5 perces frissítéssel nyomon követhetık a mőholdak. A pályára állítástól (2002. március 17. 12 óra) pedig folyamatosan számolja, mennyi ideje keringenek a mőholdak. 2005. október 21-én 22 órakor 1314 napja és 10 órája keringtek!
33
