Matematika (a fyzika) schovan´a za GPS Michal Bulant Masarykova univerzita Pˇr´ırodovˇ edeck´ a fakulta ´ Ustav matematiky a statistiky
Brno, 8. bˇrezna 2012
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
1 / 24
Global Positioning system minim´alnˇe 27 satelit˚ u (24 aktivn´ıch po 4 rovnomˇernˇe rozm´ıstˇen´e na 6 orbit´aln´ıch drah´ach, 3 z´aloˇzn´ı) – http://www.nstb.tc.faa.gov/Full_WaasSatelliteStatus.htm v´yˇska cca 19 300 km na povrchem Zemˇe, cca 2 obˇehy dennˇe z kaˇzd´eho m´ısta na Zemi viditeln´ych 4–12 satelit˚ u od 1. kvˇetna 2000 zruˇseno umˇel´e zkreslov´an´ı dat (SA – selective availability)
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
2 / 24
V´ypoˇcet pozice – u´vod Satelity ob´ıhaj´ıc´ı (nejde o stacion´arn´ı druˇzice) Zemi vys´ılaj´ı zpr´avy obsahuj´ıc´ı: ˇcas vysl´an´ı zpr´avy, polohu satelitu, syst´emovou informaci o stavu a (pˇribliˇzn´e) pozici ostatn´ıch satelit˚ u. Z tˇechto informac´ı chce pˇr´ıjemce (GPS pˇrij´ımaˇc) odvodit informaci o sv´e poloze.
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
3 / 24
V´ypoˇcet pozice Pˇrij´ımaˇc na z´akladˇe polohov´e a ˇcasov´e informace [xi , yi , zi , ti ] od alespoˇ n 3(4) satelit˚ u vypoˇcte svoji zd´anlivou vzd´alenost ri od jednotliv´ych vys´ılaˇc˚ u (pseudorange) za pˇredpokladu, ˇze se sign´al ˇs´ıˇr´ı rychlost´ı svˇetla (odhadnˇete, jak dlouho let´ı sign´al). Vypoˇcten´a vzd´alenost od satelitu spolu s jeho polohou pˇri vysl´an´ı sign´alu ud´av´a sf´eru (povrch koule), na n´ıˇz pˇrij´ımaˇc leˇz´ı. Pr˚ useˇc´ıkem takov´ych dvou sf´er je pak kruˇznice, obsahuj´ıc´ı dan´y bod.
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
4 / 24
V´ypoˇcet pozice – pokraˇcov´an´ı
Pr˚ useˇc´ıkem tˇret´ı sf´ery s touto kruˇznic´ı jsou pak (obvykle) 2 body. V´yslednou pozici je pak moˇzn´e urˇcit jako: ten z pr˚ useˇc´ık˚ u, kter´y je bl´ıˇze povrchu Zemˇe (v obvykl´em pˇr´ıpadˇe GPS pˇrij´ımaˇce v autˇe ˇci v ruce) ten z pr˚ useˇc´ık˚ u, kter´y je bl´ıˇze ˇ ctvrt´ e sf´ eˇre – v tomto pˇr´ıpadˇe je rovnˇeˇz moˇzn´e pomoc´ı GPS urˇcit nadmoˇrskou v´yˇsku, v n´ıˇz se pˇrij´ımaˇc pohybuje.
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
5 / 24
Koneˇcnˇe sl´ıben´a matematika
Pro zjednoduˇsen´ı v´ypoˇct˚ u je moˇzn´e bez u ´jmy na obecnosti zvolit kart´ezskou soustavu souˇradnic tak, ˇze stˇredy sf´er (tj. pozice vys´ılaj´ıc´ıch satelit˚ u) jsou v rovinˇe xy (tj. z = 0), jeden ze stˇred˚ u d´ale um´ıst´ıme v poˇc´atku a druh´y na ose x. Uvaˇzujme tedy tˇri sf´ery se stˇredy v bodech [0, 0, 0], [u, 0, 0], [v , w , 0] a polomˇery r1 , r2 , r3 a dostaneme tak pro hledanou pozici [x, y , z] rovnice x 2 + y 2 + z 2 = r12 (x − u)2 + y 2 + z 2 = r22 (x − v )2 + (y − w )2 + z 2 = r32
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
6 / 24
Koneˇcnˇe sl´ıben´a matematika x 2 + y 2 + z 2 = r12 (x − u)2 + y 2 + z 2 = r22 (x − v )2 + (y − w )2 + z 2 = r32 Zkuste si soustavu vyˇreˇsit! Odeˇcten´ım 2. rovnice od prvn´ı a snadnou 1 u ´pravou dostaneme x = 2u (r12 − r22 + u 2 ), odkud po dosazen´ı za x do prvn´ı rovnice dostaneme vztah r12 −
(r12 − r22 + u 2 )2 = y 2 + z 2. 4u 2
Podm´ınkou pro ˇreˇsitelnost (tj. pro to, ˇze se prvn´ı dvˇe sf´ery v˚ ubec prot´ınaj´ı) 2 2 2 2 2 je 2ur1 ≥ r1 − r2 + u , neboli r2 ≥ (u − r1 ) , ˇci r1 + r2 ≥ u ≥ r1 − r2 (tuto podm´ınku lze samozˇrejmˇe takˇrka ihned vidˇet z obr´azku). Pˇri splnˇen´ı odvozen´e podm´ınky jiˇz vypoˇcteme i souˇradnici y pomoc´ıq dosazen´ı do tˇret´ı rovnice. Souˇradnici z pak lze dopoˇc´ıtat napˇr. jako z = ± r12 − x 2 − y 2 . Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
7 / 24
Jak ale poˇc´ıtat prakticky odmocniny? V d˚ usledku je tˇreba ˇreˇsit neline´arn´ı soustavu rovnic o v´ıce nezn´am´ych – jiˇz jsme uk´azali jeden zp˚ usob, jak´ym ji lze pˇrev´est na postupn´e ˇreˇsen´ı rovnic o jedn´e nezn´am´e. Newton-Raphsonova metoda je iterativn´ı metoda na hled´an´ı koˇren˚ u re´aln´ych funkc´ı (obecnˇe v´ıce promˇenn´ych).
Newtonova metoda S touto metodou pˇriˇsel Newton kolem roku 1670 a vysvˇetlil ji na pˇr´ıkladu rovnice x 3 − 2x − 5 = 0. Jeden z koˇren˚ u je bl´ızko 2, poloˇzil tedy x = 2 + p a dosazen´ım do rovnice dostal vztah pro p: p 3 + 6p 2 + 10p − 1 = 0. 1 Protoˇze je ale p mal´e, je moˇzn´e zanedbat ˇcleny p 3 , 6p 2 , odkud p = 10 . To samozˇrejmˇe nen´ı pˇresn´e ˇreˇsen´ı, jde ale o dalˇs´ı zpˇresnˇen´ı, m˚ uˇzeme nyn´ı ps´at x = 2,1 + q, dostat tak dalˇs´ı aproximaci x = 2,0946 atd. Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
8 / 24
Jak ale poˇc´ıtat prakticky odmocniny? Ukaˇzme zde pro ilustraci pouˇzit´ı t´eto metody pro odvozen´ı elegantn´ıho postupu v´ypoˇctu druh´e odmocniny (tento postup je zn´am jako Babyl´onsk´a metoda ˇci jako Heron˚ uv vzorec1 ). 1
Mˇejme d´anu diferencovatelnou funkci f (x) a aproximaci jej´ıho koˇrene x0 .
2
Postupnˇe poˇc´ıtejme dalˇs´ı iterace pomoc´ı vztahu xn+1 = xn −
Pro v´ypoˇcet druh´e odmocniny z a (tj. hled´an´ı koˇrene funkce f (x) = x 2 − a) tak dost´av´ame iteraˇcn´ı postup xn+1 = 21 (xn +
f (xn ) f 0 (xn ) .
a xn ).
Tato metoda se d´a analogicky pouˇz´ıt pˇri optimalizaci, kde m´ısto koˇrene hled´ame ˇreˇsen´ı rovnice f 0 (x) = 0. 1 To samozˇrejmˇe neznamen´ a, ˇze Newton mˇel nˇeco spoleˇcn´eho s d´ avn´ymi Babyl´ on ˇany, jeho metoda je obecnˇejˇs´ı. Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
9 / 24
Pˇr´ıklad
√ Vypoˇctˇeme 12 s x0 = 3: x1 = √ pˇritom 12 ≈ 3,46410.
3+4 2 ,
x2 =
7/2+24/7 2
= 97/28 ≈ 3,46429,
Anal´yza efektivity Newtonovy metody Pomoc´ı Taylorovy vˇety lze v nˇejak´em okol´ı xn ps´at f (x) = f (xn ) + f 0 (xn )(x − xn ) +
1 00 f (α)(x − xn )2 , 2!
nuj´ıc´ı f (x) = 0, lze po kde α je mezi xn a x. Protoˇze hled´ame x splˇ vydˇelen´ı f 0 (xn ) vztah upravit na f (xn ) f 00 (α) + (x − x ) = − (x − xn )2 , n f 0 (xn ) 2f 0 (xn ) a tedy x − xn+1 = − Michal Bulant (PˇrF MU)
f 00 (α) (x − xn )2 . 2f 0 (xn )
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
10 / 24
Newtonova metoda – pˇr´ıklad, kdy nefunguje
Pˇr´ıklad Pˇr´ıkladem funkce, jej´ıˇz koˇren tato metoda nenajde, ani kdyˇz zaˇcneme √ sebebl´ıˇze, je f (x) = 3 x. Zde totiˇz dostaneme 1/3
xn+1 = xn+1 = xn −
Michal Bulant (PˇrF MU)
f (xn ) xn = xn − −2/3 = −2xn . 0 1 f (xn ) 3 xn
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
11 / 24
Zobecnˇen´ı na pˇr´ıpad v´ıce promˇenn´ych Zobecnˇen´ı na (napˇr.) k rovnic o k nezn´am´ych je relativnˇe pˇr´ımoˇcar´e: xn+1 = xn − JF (xn )−1 · F (xn ), kde JF je Jacobi´an zobrazen´ı F . V´ypoˇcet jeho inverze je ale ˇcasovˇe velmi n´aroˇcn´a operace, proto se ˇcasto m´ısto toho vyuˇz´ıv´a ˇreˇsen´ı pˇr´ısluˇsn´e soustavy line´arn´ıch rovnic, v´ypoˇcet zobecnˇen´e inverze, pˇri v´ıce neˇz k rovnic´ıch metoda nejmenˇs´ıch ˇctverc˚ u metoda sdruˇzen´ych gradient˚ u pro ˇreˇsen´ı pˇr´ısluˇsn´e soustavy, r˚ uzn´ych tzv. kvazi-newtonovsk´ych metod, vyuˇz´ıvaj´ıc´ıch pouze pˇribliˇzn´eho Hessi´anu (napˇr. BFGS) – viz napˇr. http://demonstrations.wolfram.com/ MinimizingTheRosenbrockFunction/. Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
12 / 24
Fyzika a praxe n´am to trochu zkomplikuje
Do ide´aln´ıho stavu uk´azan´eho dˇr´ıve se n´am ale vloud´ı v´ıce ˇci m´enˇe z´avaˇzn´e chyby:
2
Satelity disponuj´ı vysoce pˇresn´ymi atomov´ymi hodinami, to ale naˇse kapesn´ı GPSka neum´ı (st´ala by ˇr´adovˇe mili´ ony). ˇıˇr´ı se sign´al skuteˇcnˇe rychlost´ı svˇetla i pˇri pr˚ S´ uchodu ionosf´erou?
3
Sign´al se odr´aˇz´ı od r˚ uzn´ych ter´enn´ıch pˇrek´aˇzek, budov apod.
4
Do hry velmi z´asadnˇe vstupuje i speci´aln´ı a obecn´a teorie relativity.
1
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
13 / 24
Zdroje chyb GPS Error Source Ionosphere Troposphere Satellite Clock Synchronization Electronic Noise Multipath Error Satellite Position (Ephemeris) Intentional Degradation Net RMS error Typical Geometric Error (GDOP) Final RMS error (Net x GDOP) Actual Typical Error
Typical or Maximum Error 10 Meters 1 Meter 1 Meter 2 Meters 0.5 Meters 1 Meter 0 Meters 10 Meters 4 40 meters 10 meters
Zdroj: http://www.pdhcenter.com/courses/l116/l116content.htm
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
14 / 24
Jak se vyrovnat s chybami – hodiny v pˇrij´ımaˇci S nepˇresnost´ı levn´ych hodin v GPS pˇrij´ımaˇci se vyrovn´ame pomˇernˇe snadno – k tomu n´am slouˇz´ı pr´avˇe ˇctvrt´y (a pˇr´ıpadnˇe dalˇs´ı) satelit, kter´y jsme dosud ve v´ypoˇctech nepouˇzili. V praxi tak dost´av´ame ˇctyˇri nebo v´ıce rovnic o ˇctyˇrech nezn´am´ych (x, y , z, error ). Na obr´azku je pro zjednoduˇsen´ı uk´az´an 2D pˇr´ıpad, kde hodiny v pˇrij´ımaˇci jsou zpoˇzdˇeny o 0,5 s.
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
15 / 24
Jak se vyrovnat s chybami – hodiny v pˇrij´ımaˇci Pokud je vidˇet v´ıce neˇz ˇctyˇri satelity, m´ame tzv. pˇreurˇcen´y syst´em rovnic a do hry vstupuje moˇznost vybrat si z nˇekolika moˇznost´ı tu nejlepˇs´ı – v takov´em pˇr´ıpadˇe se poloha aproximuje pomoc´ı metody nejmenˇs´ıch ˇctverc˚ u. Metoda slouˇz´ı k rekonstrukci funkce f z hodnot f0 , . . . , fn namˇeˇren´ych v uzlov´ych bodech a0 , . . . , an . Tuto rekonstrukci hled´ame vzhledem k dan´emu modelu – dan´e posloupnosti funkc´ı (obecnˇe v´ıce promˇenn´ych) g0 (x), . . . , gm (x), . . . – ve tvaru ym (x) =
m X
cj gj (x).
j=0
C´ılem je pˇri tom minimalizovat souˇcet ˇctverc˚ u n X
2 fi − ym (ai ) .
i=0 Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
16 / 24
Aproximace metodou nejmenˇs´ıch ˇctverc˚ u
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
17 / 24
Ukaˇzme si pouˇzit´ı t´eto metody v nejjednoduˇsˇs´ım pˇr´ıpadˇe, kdy m´ame d´ano n bod˚ u ([x1 , y1 ], . . . , [xn , yn ]) a hled´ame pˇr´ımku, kter´a nejl´epe vystihuje rozloˇzen´ı tˇechto bod˚ u. Hled´ame tedy funkci tvaru f (x) = a · x + b s nezn´am´ymi a, b ∈ R tak, aby hodnota n X (f (xi ) − yi )2 i=1
byla minim´aln´ı. S vyuˇzit´ım diferenci´aln´ıho poˇctu lze snadno odvodit n´asleduj´ıc´ı tvrzen´ı.
Vˇeta Mezi pˇr´ımkami tvaru f (x) = a · x + b m´a nejmenˇs´ı souˇcet ˇctverc˚ u vzd´alenost´ı funkˇcn´ıch hodnot v bodech x1 , . . . , xn od hodnot yi funkce splˇ nuj´ıc´ı X X X a xi2 + b xi = xi yi X X a xi + b · n = yi Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
18 / 24
Metoda nejmenˇs´ıch ˇctverc˚ u – pˇr´ıklad Pˇr´ıklad Metodou nejmenˇs´ıch ˇctverc˚ u urˇcete regresn´ı pˇr´ımku odpov´ıdaj´ıc´ı x 1 2 3 4 namˇeˇren´ym dat˚ um: y 1,5 1,6 2,1 3,0
ˇ sen´ı Reˇ Data je vhodn´e seˇradit v tabulce podle sch´ematu: x 1 2 3 4 10
y 1,5 1,6 2,1 3 8,2
xy 1,5 3,2 6,3 12 23
x2 1 4 9 16 30
Odtud 30a + 10b = 23, 10a + 4b = 8,2, a tedy a = 0,5, b = 0,8. Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
19 / 24
Jak se vyrovnat s chybami – teorie relativity GPS ukazuje jeden z nejpraktiˇctˇejˇs´ıch d˚ usledk˚ u teorie relativity – pokud bychom ji nevzali v potaz, bude metoda GPS prakticky nepouˇziteln´a. Atomov´e hodiny pracuj´ı s pˇresnost´ı na nanosekundy (ns = 10−9 s), abychom byli schopni zaruˇcit pˇresnost zjiˇstˇen´ı pozice na cca 10 m, je tˇreba umˇet urˇcit pˇresnost ˇcasu vys´ılaˇce s pˇresnost´ı cca 30 ns. Pˇritom se satelity vzhledem k Zemi pohybuj´ı rychlost´ı cca 14 000 km/h. Do hry tak vstupuje speci´aln´ı teorie relativity, nebot’ pˇrij´ımaˇc a vys´ılaˇc jsou v˚ uˇci sobˇe v pohybu, doch´az´ı ke zpomalen´ı hodin vys´ılaˇce oproti v2 42 −10 , tj. asi pozorovateli (dilatace ˇcasu) o 2c 2 ≈ 2·(3·105 )2 ≈ 10 o 7,7 µs/den. Dalˇs´ı jeˇstˇe v´yznamnˇejˇs´ı efekt pˇredstavuje obecn´a teorie relativity, kter´a implikuje, ˇze hodiny pobl´ıˇz masivn´ıho objektu (Zemˇe) jdou pomaleji neˇz hodiny vzd´alenˇejˇs´ı (d´ıky vˇetˇs´ımu zakˇriven´ı prostoroˇcasu). Z povrchu Zemˇe vid´ıme tedy satelitn´ı hodiny jdouc´ı rychleji neˇz tyt´eˇz hodiny um´ıstˇen´e na Zemi o cca 45µs za den. Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
20 / 24
Jak se vyrovnat s chybami – teorie relativity
Nezapoˇc´ıt´an´ım teorie relativity bychom tak dostali chybu v ˇr´adu 38µs za den, coˇz v d˚ usledku znamen´a cca 10km chybu v urˇcen´ı pozice. Tato chyba je opravena umˇel´ym zpomalen´ım atomov´ych hodin um´ıstˇen´ych v satelitech oproti hodin´am na Zemi (10,22999999543 MHz oproti 10,23 MHz).
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
21 / 24
Diferenci´aln´ı GPS
Jedno z mnoha moˇzn´ych vylepˇsen´ı je zaloˇzeno na myˇslence, ˇze relativnˇe bl´ızk´e pˇrij´ımaˇce podl´ehaj´ı analogick´ych atmosf´erick´ym chyb´am. D´ıky pevn´ym stanic´ım (http://www.ndblist.info/datamodes/worldDGPSdatabase.pdf), u nichˇz je s vysokou pˇresnost´ı zn´ama poloha a kter´e vys´ılaj´ı rozd´ıl mezi touto polohou a polohou vypoˇctenou na z´akladˇe informac´ı ze satelit˚ u, je moˇzn´e u ˇspiˇcov´ych DGPS pˇr´ıstroj˚ u dos´ahnout pˇresnosti v ˇr´adu centimetr˚ u.
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
22 / 24
Pˇr´ıklad na z´avˇer Pˇr´ıklad V tabulce jsou uvedena skuteˇcn´a data z nˇekolika satelit˚ u – geocentrick´e souˇradnice jsou uvedeny v metrech, ˇcas pˇrenosu sign´alu v nanosekund´ach. Vaˇsim u ´kolem je s vyuˇzit´ım vhodn´eho SW (napˇr OpenOffice Calc) urˇcit: 1
geocentrick´e souˇradnice m´ısta pozorovatele,
2
popsat skuteˇcn´e m´ısto na Zemi, kde se pozorovatel nach´azel.
ˇ sat. C. 1 2 3 4 5 6
x [m] 14177553.47 15097199.81 23460342.33 -8206488.95 1399988.07 6995655.48
Michal Bulant (PˇrF MU)
y [m] -18814768.09 -4636088.67 -9433518.58 -18217989.14 -17563734.90 -23537808.26
z [m] 12243866.38 21326706.55 8174941.25 17605231.99 19705591.18 -9927906.48
Matematika (a fyzika) schovan´ a za GPS
dt [ns] 70446329.64 75142197.81 78968497.2 69887173.01 67231182.38 80796265.09 Brno, 8. bˇrezna 2012
23 / 24
Pouˇzit´a literatura
Wikipedia, The Free Encyclopedia, www.wikipedia.org. Neil Ashby, Relativity and the Global Positioning System. Physics Today, May 2002.
Dˇekuji za pozornost!
Michal Bulant (PˇrF MU)
Matematika (a fyzika) schovan´ a za GPS
Brno, 8. bˇrezna 2012
24 / 24
