ČESKÉ VYSOKÉ UČENÍ TECHNICKÉ V PRAZE FAKULTA STAVEBNÍ – OBOR GEODÉZIE A KARTOGRAFIE KATEDRA VYŠŠÍ GEODÉZIE název předmětu
Vyšší geodézie 1 úloha/zadání
3/3 školní rok
2010/11
název úlohy
GPS - výpočet polohy stanice pomocí dvojitých diferencí semestr
1
skupina
NG1-88
zpracoval
datum
Jan Dolista
22. 11.
klasifikace
GPS - výpočet polohy stanice pomocí dvojitých diferencí Zadání: Máte k dispozici observační a navigační soubory formátu RINEX, pořízené pro body 4001, 4003 a 4004. V první části úlohy vypočítejte polohu a opravu hodin přijímače referenčního bodu 4004 na základě časových epoch 𝑡1 , 𝑡2 , 𝑡3 observačního souboru RINEX. Výpočet proveďte zpracováním kódových měření, pracujte s pseudovzdálenostmi získanými v P2-kódu. Polohy družic a opravy družicových hodin určujte na základě příslušného navigačního souboru formátu RINEX (pokud nejsou v daném souboru k dispozici nejbližší dřívější efemeridy - nejbližší nižší čas efemerid, tak použijte nejbližší pozdější). Jako přibližnou polohu bodu 4004 vstupující do vyrovnání použijte přibližné souřadnice přijímače udané v observačním souboru formátu RINEX. V druhé části úlohy vypočítejte polohu připojovaného bodu 4001/4003 relativně k bodu 4004. Pracujte opět s pseudovzdálenostmi získanými v P2-kódu ve stejných časových epochách a výpočet proveďte zpracováním dvojitě diferencovaných měření. Při výpočtu berte v úvahu družice společné pro referenční a připojovaný bod. Jako rychlost světla použijte hodnotu 𝑐 = 299792458.0𝑚𝑠−1 , jako geocentrickou gravitační konstantu 𝐺𝑀 = 398, 6005 · 1012 𝑚3 𝑠−2 a jako úhlovou rychlost rotace Země 𝜔𝑒 = 7292115, 1467 · 10−11 𝑠−1 .
Číselné zadání 3: číslo zadání referenční bod 3 4004
připojovaný bod 4003
epocha 𝑡1 [s] 48150
epocha 𝑡2 [s] 48750
epocha 𝑡3 [s] 49350
Epochy 𝑡1 , 𝑡2 , 𝑡3 jsou zadány v sekundách od začátku dne.
Vypracování: Veškeré výpočty byly provedeny v programu Octave.
1 1.1
Poloha bodu 4004 Načtení měřené pseudovzdálenosti a přibližných souřadnic stanice
Jako první byly zadané časy epoch převedeny na hodiny, minuty a vteřiny v rámci daného dne. Pro tyto časy byla z observačního RINEXu vybrána měření na všechny družice NAVSTAR GPS. Pro další zpracování bylo načteno PRN družice a měřená pseudovzdálenost získaná z P2-kódu. Z hlavičky souboru byly převzaty přibližné souřadnice určovaného bodu [𝑋0 , 𝑌0 , 𝑍0 ]. Pro tento účel byla v programu Octave vytvořena funkce, která pro zadaný observační soubor a čas epochy vrací vektory výše zmíněných hodnot (PRN, pseudovzdálenost a přibližné souřadnice). Funkce ještě vyžaduje pomocný parametr a to rok měření.
1.2
Poloha družice v době vyslání signálu
Na základě již načtených pseudovzdáleností byl pro každou z družic určen přibližný čas vyslání signálu: 𝑅 𝑡𝑖1 = 𝑡𝑒𝑝𝑜𝑐ℎ𝑦 − , 𝑐 kde 𝑡𝑒𝑝𝑜𝑐ℎ𝑦 je zadaný čas měření, 𝑅 měřená pseudovzdálenost a 𝑐 rychlost světla. Pro tento čas a každou z družic byl v navigačním RINEXu vyhledán záznam s nejbližším nižším časem, pokud takový záznam nebyl nalezen, byl použit nejbližší vyšší čas měření. Pro vyhledání záznamu a výpočet polohy byla opět vytvořena funkce v programu Octave, jejíž vstupní parametry
jsou PRN, přibližný čas vyslání signálu, soubor navigačního RINEXu a pomocný parametr rok měření. Funkce po nalezení příslušného záznamu ve tvaru: PRN — 𝐶𝑢𝑐 𝑇 𝑜𝑒 𝐼0 𝐼˙ .. .
ROK 𝐶𝑟𝑠 𝑒 𝐶𝑖𝑐 𝐶𝑟𝑐 —
MĚSÍC ∆𝑛 𝐶𝑢𝑐 𝑙0 𝜔0 —
DEN 𝑀0 √ 𝑎 𝐶𝑖𝑠 ˙ Ω —
HODINA
MINUTA
SEKUNDA
𝛿𝑐
𝜕𝛿𝑐 𝜕𝑡
𝜕 2 𝛿𝑐 𝜕𝑡2
načte parametry dráhy družice do proměnných a spočte souřadnice družice v systému WGS-84 a chybu hodin. Tyto hodnoty pak vrátí jako výstup. 1.2.1
Přesný čas vyslání signálu a chyba hodin družice
Jako první je pro danou družici určena chyba hodin v přibližném čase vyslání signálu 𝑡𝑖1 : 𝛿𝑠 = 𝛿𝑐 +
𝜕𝛿𝑐 𝜕 2 𝛿𝑐 · (𝑡𝑖1 − 𝑡0 ) + · (𝑡𝑖1 − 𝑡0 )2 𝜕𝑡 𝜕𝑡2
Dále je zpřesněn čas vyslání signálu: 𝑡1 = 𝑡𝑖1 − 𝛿𝑠 Vzhledem k malým rozdílům časů 𝑡1 a 𝑡𝑖1 již nebyla dále zpřesňována chyba hodin družice. Dalším krokem byl výpočet Keplerovských oskulačních elementů: 1.2.2
Střední anomálie 𝑀 = 𝑀0 + 𝑛(𝑡1 − 𝑡0 ) + ∆𝑛(𝑡1 − 𝑡0 )
1.2.3
Excentrická anomálie 𝐸 = 𝑀 + 𝑒 sin 𝐸
Excentrická anomálie byla určena iteračně, kdy v první iteraci je hodnota excentrické anomálie volena 𝐸0 = 𝑀 + 𝑒 sin 𝑀 . V dalších iteracích je hodnota anomálie dána vztahem 𝐸𝑖 = 𝑀 + 𝑒 sin 𝐸𝑖−1 . Výpočet je opakován, dokud rozdíl ve dvou po sobě jdoucích iteracích není menší než dvojnásobek přesnosti výpočtu (dáno možnostmi počítače). 1.2.4
Pravá anomálie
(︃ √
1+𝑒 𝐸 𝑣 = 2𝑎𝑡𝑎𝑛 √ 𝑡𝑔 1−𝑒 2 1.2.5
)︃
Argument šířky 𝑢0 = 𝜔0 + 𝑣 𝜔 = 𝜔0 + 𝐶𝑢𝑐 cos (2𝑢0 ) + 𝐶𝑢𝑠 sin (2𝑢0 ) 𝑢=𝜔+𝑣
1.2.6
Průvodič 𝑟0 = 𝑎(1 − 𝑒 cos 𝐸) 𝑟 = 𝑟0 + 𝐶𝑟𝑐 cos (2𝑢) + 𝐶𝑟𝑠 sin (2𝑢)
Po výpočtu argumentu šířky a délky průvodiče byly určeny souřadnice družice v systému, kde oběžná dráha družice je v rovině xy, a osa x směřuje do výstupního uzlu (průsečík dráhy s rovinou rovníku). ⎞ ⎛ cos 𝑢 ⎟ ⎜ XS = 𝑟 ⎝ sin 𝑢 ⎠ 0 1.2.7
Sklon ˙ 1 − 𝑡0 ) 𝐼 = 𝐼0 + 𝐶𝑖𝑐 cos (2𝑢) + 𝐶𝑖𝑠 sin (2𝑢) + 𝐼(𝑡
1.2.8
Délka výstupního uzlu ˙ − 𝜔𝐸 )(𝑡1 − 𝑡0 ) − 𝑇 𝑜𝑒𝜔𝐸 𝑙 = 𝑙0 + (Ω
1.3
Rotace do WGS
Posledním krokem výpočtu souřadnic družice je rotace do systému WGS. Ta probíhá ve dvou krocích. Nejprve sklopením kolem osy x o úhel 𝐼 v matematicky záporném smyslu. Tím je ztotožněna osa z s osou Z systému WGS. Druhým krokem je pak rotace kolem osy Z o úhel 𝑙 v matematicky záporném smyslu. XWGS = RZ (−𝑙)RX (−𝐼)XS , kde ⎛
⎞
cos (−𝑙) sin (−𝑙) 0 ⎜ ⎟ RZ (−𝑙) = ⎝ − sin (−𝑙) cos (−𝑙) 0 ⎠ 0 0 1
1.4
⎛
⎞
1 0 0 ⎜ ⎟ RX (−𝐼) = ⎝ 0 cos (−𝐼) sin (−𝐼) ⎠ 0 − sin (−𝐼) cos (−𝐼)
Určení polohy bodu vyrovnáním MNČ
Pomocí obou funkcí byly získany přibližné souřadnice stanice 𝑋𝑅0 , 𝑌𝑅0 , 𝑍𝑅0 a měřené pseudovzdá𝑆 na všechny družice dostupné v daných etapách. Ke každé pseudovzdálenosti byly také lenosti 𝑅𝑅 určeny souřadnice družice v době vyslání signálu 𝑋 𝑆 , 𝑌 𝑆 , 𝑍 𝑆 a chyba hodin družice 𝛿𝑠 v tomto okamžiku. Pro výpočet vyrovnání byla vytvořena další funkce v programu Octeve, jejímž vstupem je vektor přibližných souřadnic stanice, matice, která na každém řádku obsahuje: 𝑆 𝑋𝑆 𝑌 𝑆 𝑍𝑆 𝛿 , PRN 𝑅𝑅 𝑠 posledním vstupem je pak vektor obsahující počty měření v jednotlivých etapách. 1.4.1
Přibližná geometrická vzdálenost mezi určovaným bodem a družicí
Z těchto hodnot lze určit přibližnou geometrickou vzdálenost mezi družicí a stanicí: 𝜌𝑆𝑅0 = 1.4.2
√︁
(𝑋𝑅0 − 𝑋 𝑆 )2 + (𝑌𝑅0 − 𝑌 𝑆 )2 + (𝑍𝑅0 − 𝑍 𝑆 )2
Vektor redukovaných měření
Vektor l obsahuje redukovaná měření pro všechny použité družice ve všech třech epochách. Jednotlivé prvky vektoru l jsou určeny 𝑆 𝑙𝑖 = 𝜌𝑆𝑅0 − 𝑅𝑅 − 𝑐𝛿 𝑆 , 𝑆 je měřená pseudovzdálenost a 𝛿 𝑆 chyba hodin kde 𝜌𝑆𝑅0 je přibližná geometrická vzdálenost, 𝑅𝑅 družice.
1.4.3
Derivace funkčního vztahu podle neznámých
Měřená pseudovzdálenost je vzhledem k neznámým souřadnicím a chybě hodin přijímače vyjádřena vztahem: 𝜌𝑆𝑅 =
√︁
(𝑋𝑅 − 𝑋 𝑆 )2 + (𝑌𝑅 − 𝑌 𝑆 )2 + (𝑍𝑅 − 𝑍 𝑆 )2 + 𝑐𝛿𝑅 − 𝑐𝛿 𝑆
Derivace tohoto vztahu podle souřadnic: 𝑎𝑆𝑅0 =
𝑋𝑅0 − 𝑋 𝑆 𝜌𝑆𝑅0
𝑏𝑆𝑅0 =
𝑌𝑅0 − 𝑌 𝑆 𝜌𝑆𝑅0
𝑐𝑆𝑅0 =
𝑍𝑅0 − 𝑍 𝑆 𝜌𝑆𝑅0
Derivace tohoto vztahu podle chyby hodin přijímače: 𝑑𝑆𝑅0 = 1 Derivováno bylo podle chyby hodin přijímače násobené rychlostí světla, tento fakt je nutné v dalším výpočtu zohlednit. 1.4.4
Matice plánu
submatice pro 1. epochu: 𝑎1𝑅0 𝑏1𝑅0 𝑐1𝑅0 1 0 0 ⎜ .. .. .. .. .. ⎟ , A1 = ⎝ ... . . . . . ⎠ 𝑎𝑘𝑅0 𝑏𝑘𝑅0 𝑐𝑘𝑅0 1 0 0 ⎞
⎛
kde 𝑘 je počet družic v 1. epoše submatice pro 2. epochu: 𝑎1𝑅0 𝑏1𝑅0 𝑐1𝑅0 0 1 0 ⎜ .. .. .. .. .. .. ⎟ , A2 = ⎝ . . . . . . ⎠ 𝑝 𝑝 𝑝 𝑎𝑅0 𝑏𝑅0 𝑐𝑅0 0 1 0 ⎞
⎛
kde 𝑝 je počet družic ve 2. epoše submatice pro 3. epochu: 𝑎1𝑅0 𝑏1𝑅0 𝑐1𝑅0 0 0 1 ⎜ .. .. .. .. .. ⎟ , A3 = ⎝ ... . . . . . ⎠ 𝑞 𝑞 𝑞 𝑎𝑅0 𝑏𝑅0 𝑐𝑅0 0 0 1 ⎞
⎛
kde 𝑞 je počet družic ve 3. epoše výsledná matice plánu: ⎛
⎞
𝐴1 ⎜ ⎟ A = ⎝ 𝐴2 ⎠ 𝐴3 1.4.5
Opravy přibližných hodnot neznámých (vlastní vyrovnání)
Opravy přibližných hodnot souřadnic a chyb hodin přijímače byly vypočteny vyrovnáním MNČ. (︁
dh = − A𝑇 A
)︁−1
· A𝑇 l
1.4.6
Vyrovnané hodnoty souřadnic bodu a chyb hodin přijímače
Vyrovnané hodnoty neznámých se vypočtou přičtením oprav k přibližným hodnotám. Přibližné hodnoty chyb hodin přijímače jsou 𝛿𝑅1 = 𝛿𝑅2 = 𝛿𝑅3 = 0. ⎛ ⎜ ⎜ ⎜ ⎜ h=⎜ ⎜ ⎜ ⎝
𝑋𝑅0 𝑌𝑅0 𝑍𝑅0 0 0 0
⎞ ⎟ ⎟ ⎟ ⎟ ⎟ + dh ⎟ ⎟ ⎠
Jelikož byla v matici plánu použita derivace podle chyb hodin přijímače násobené rychlostí světla, (𝑐𝛿𝑅 ) bylo nutné opravy chyb hodin přijímače vydělit konstantou 𝑐. h= 1.4.7
(︁
𝑋𝑅 𝑌𝑅 𝑍𝑅 𝛿𝑅1 𝛿𝑅2 𝛿𝑅2
)︁𝑇
Vektor vyrovnaných měření
S využitím vyrovnaných hodnot neznámých byl vypočten vektor vyrovnaných měření lvyr : 𝑙𝑖𝑣𝑦𝑟 1.4.8
=
√︁
(𝑋𝑅 − 𝑋 𝑆 )2 + (𝑌𝑅 − 𝑌 𝑆 )2 + (𝑍𝑅 − 𝑍 𝑆 )2 + 𝑐𝛿𝑅 − 𝑐𝛿 𝑆
I. a II. výpočet oprav
Opravy byly vypočteny dvakrát. Jednou v rámci vyrovnání a podruhé jako rozdíl pseudovzdáleností vypočtených z vyrovnaných neznámých a měřených pseudovzdáleností. vI = A · dh + l vII = lvyr − lP2 , kde lP2 je vektor obsahující měřené pseudovzdálenosti. 1.4.9
Porovnání I. a II. výpočtu oprav
I. a II. výpočet oprav byly vzájemně porovnány. Jelikož rozdíl řádově 10−6 byl považován za dostatečně malý, nebyla provedena další iterace. vI − vII
1.5
Charakteristiky přesnosti
Na závěr byly vypočteny charakteristiky přesnosti. 1.5.1
Aposteriorní střední jednotková chyba √︃
𝑚0 =
vT v , 𝑛−𝑚
kde n=k+p+q ... počet družic, m=6 ... počet neznámých. K výpočtu byl použit I. výpočet oprav. 1.5.2
Kovarianční matice neznámých (︁
Qdh = AT A
)︁−1
1.5.3
Střední chyby neznámých mdh = 𝑚0 ·
√︁
𝑑𝑖𝑎𝑔(Qdh )
Střední chyby chyb hodin přijímače jsou opět vzhledem k použité derivaci v matici plánu násobeny rychlostí světla 𝑐.
1.6
Číselné výsledky pro bod 4004
Přibližné souřadnice určovaného bodu: 𝑋0 = 3970656.0𝑚 𝑌0 = 1018570.3𝑚 𝑍0 = 4870425.3𝑚 Střední aposteriorní jednotková chyba: 𝑚0 = 4.1 Vyrovnané hodnoty souřadnic a chyb hodin přijímače a jejich střední chyby: 𝑋 = 3970663.1𝑚 𝑚𝑋 = 2.6𝑚 𝑌 = 1018585.9𝑚 𝑚𝑌 = 1.7𝑚 𝑍 = 4870434.1𝑚 𝑚𝑍 = 2.8𝑚 𝛿𝑅1 = 8.0393643𝑒 − 07𝑠 𝑚𝛿𝑅1 = 9.170𝑒 − 09𝑠 𝛿𝑅2 = −0.00022616521𝑠 𝑚𝛿𝑅2 = 9.122𝑒 − 09𝑠 𝛿𝑅3 = −0.00045191807𝑠 𝑚𝛿𝑅3 = 9.027𝑒 − 09𝑠
2
Poloha bodu 4003
Pro výpočet souřadnic bodu 4003 pomocí druhých diferencí byly souřadnice bodu 4004 zafixovány na hodnotách vypočtených v první části úlohy: 𝑋4004 = 3970663.1𝑚 𝑌4004 = 1018585.9𝑚 𝑍4004 = 4870434.1𝑚
2.1
Načtení měřené pseudovzdálenosti a přibližných souřadnic stanice
Pomocí výše popsané funkce v programu Octave byla z observačních RINEXů pro stanice na bodech 4004 a 4003 načtena měření pro dané časy epoch 𝑡1 , 𝑡2 , 𝑡3 . Použita byla opět kódová měření (P-kód) na nosné vlně L2. Současně byly načteny přibližné souřadnice obou stanic. Jelikož ale pro stanici 4004 jsou již známé „přesné“ souřadnice nejsou přibližné souřadnice potřeba. Souřadnice bodu 4004 jsou určeny se střední chybou 2m a tedy označení „přesné“ je zavádějící, ale pro další výpočet jsou fixovány.
2.2
Poloha družice v době vyslání signálu
Pro každé měření byl opět určen přibližný okamžik vyslání signálu: 𝑡1 = 𝑡𝑖1 − 𝛿𝑠 Druhou z již zmíněných funkcí byly pro tento okamžik vypočteny z navigačních RINEXů polohy všech družic, z nichž byl přijat signál a tedy měřena pseudovzdálenost. Vzhledem k rozdílné poloze bodů je měřená pseudovzdálenost na družici se stejným PRN různá a tedy i okamžik vyslání signálu je odlišný. Odlišná je tedy i poloha družice v okamžiku vyslání. Tyto rozdíly v poloze družice dosahují řádově decimetrů, avšak pro další výpočet byly uvažovány. Oproti tomu chod hodin družice dosahuje velmi malých hodnot. Chyba hodin družice v okamžiku vyslání pro bod 4004 se od chyby hodin družice v okamžiku vyslání pro bod 4003 liší řádově o 10−14 𝑠 a lze jí pro oba okamžiky považovat za shodnou, čehož je dále ve výpočtu využito.
Z přibližných souřadnic bodu 4003 resp. fixovaných souřadnic bodu 4004 a souřadnic družice byly vypočteny přibližné geometrické vzdálenosti: 𝜌𝑆40030 =
√︁
(𝑋40030 − 𝑋 𝑆 )2 + (𝑌40030 − 𝑌 𝑆 )2 + (𝑍40030 − 𝑍 𝑆 )2
resp. 𝜌𝑆4004 =
√︁
(𝑋4004 − 𝑋 𝑆 )2 + (𝑌4004 − 𝑌 𝑆 )2 + (𝑍4004 − 𝑍 𝑆 )2
Načtené seznamy družic v každé epoše byly porovnány a pro další výpočet byly použity pouze společné družice resp. měření na ně pro obě stanoviska.
2.3
Výpočet 1. diferencí
V dalším kroku byly vypočteny 1. diference měřených pseudovzdáleností: 𝑆 𝑆 𝑆 𝑅4003,4004 = 𝑅4004 − 𝑅4003
a zároveň i 1. diference přibližných geometrických vzdáleností: 𝜌𝑆40030 ,4004 = 𝜌𝑆4004 − 𝜌𝑆40030 , kde 𝑆 je PRN družice společné pro obě stanoviska. Při přechodu na 1. diference dojde k vyloučení chyby hodin družice, neboť je považována za přibližně stejnou pro měření z obou stanovisek, jak bylo uvedeno výše. Dále je odstraněn vliv inosférické refrakce, který lze pro obě stanoviska považovat za stejný (s vysokou přesností) a částečně je vyloučen i vliv troposférické refrakce. Troposférickou refrakci nelze pro oba přijímače považovat za stejnou, neboť největší vliv mají vlastnosti spodní vrstvy atmosféry, tedy podmínky v okolí přijímačů, které se mohou i značně lišit. Podle přibližných souřadnic jsou však oba přijímače velmi blízko sebe (několik metrů) a tedy chyba při považování vlivu troposféry za stejný nebude příliš velká.
2.4
Výpočet 2. diferencí
Poté byly vypočteny 2. diference měřených pseudovzdáleností a to odečtením měřené pseudovzdálenosti od jedné referenční (libovolně zvolené): 𝑘𝑙 𝑘 𝑙 𝑅4003,4004 = 𝑅4003,4004 − 𝑅4003,4004 .
Zároveň byly vypočteny i 2. diference přibližných geometrických vzdáleností: 𝑘 𝑙 𝜌𝑘𝑙 40030 ,4004 = 𝜌40030 ,4004 − 𝜌40030 ,4004 ,
kde 𝑘 je PRN referenční družice a 𝑙 je PRN ostatních družic v epoše. Jako referenční byla zvolena vždy první družice v dané epoše. Přechodem na 2. diference je vyloučena chyba hodin přijímače.
2.5
Matice vah
Vzhledem k tomu, že při přechodu na druhé diference je zvolena jedna družice jako referenční, je toto měření ve výpočtu použito vícekrát a měření jsou tedy závislá. Proto je nutné zavést do výpočtu matici vah. Jelikož první diference jsou nezávislé, lze matici vah 1. diferencí považovat za diagonální: ⎛ ⎜ ⎜ ⎝
PDi = ⎜
1 0 0 0
0 1 0 0
0 0 1 0
0 0 0 1
⎞ ⎟ ⎟ ⎟, ⎠
Rozměr matice odpovídá počtu 1. diferencí v rámci jedné etapy. Vztah mezi 1. a 2. diferencí lze zapsat pomocí matice F, což je matice derivací funkčního vztahu podle 1. diferencí. Obdobně jako matice plánu při vyrovnání. ⎞
⎛
1 −1 0 0 ⎟ ⎜ F = ⎝ 1 0 −1 0 ⎠ , 1 0 0 −1 Počet sloupců matice F odpovídá počtu 1. diferencí a počet řádku počtu 2. diferencí, tedy o 1 méně než 1. diferencí. Matice vah 2. diferencí je pak podle zákona přenášení vah: (︁
PDDi = F · (PDi )−1 · F𝑇
)︁−1
Druhé diference mezi epochami jsou však nezávislé a tak matici vah pro všechny epochy lze sestavit jako diagonální matici, kde na diagonále jsou submatice vah pro jednotlivé etapy a mimo diagonálu 0. ⎛ ⎞ 𝑃𝐷𝐷1 0 0 ⎜ ⎟ 𝑃𝐷𝐷2 0 ⎠, PDD = ⎝ 0 0 0 𝑃𝐷𝐷3
2.6 2.6.1
Určení polohy bodu vyrovnáním MNČ Vektor redukovaných měření
Vektor l obsahuje redukovaná měření (2. diference). Jednotlivé prvky vektoru l jsou určeny 𝑘𝑙 𝑙𝑖 = 𝜌𝑘𝑙 40030 ,4004 − 𝑅4003,4004 , 𝑘𝑙 kde 𝜌𝑘𝑙 40030 ,4004 je 2. diference přibližné geometrické vzdálenosti a 𝑅4003,4004 je 2. diference měřené pseudovzdálenosti.
2.6.2
Derivace funkčního vztahu podle neznámých
Druhá diference pseudovzdálenosti je vzhledem k neznámým souřadnicím přijímače vyjádřena vztahem: 𝑘𝑙 𝑅4003,4004 = (𝜌𝑘4004 − 𝜌𝑙4004 ) + (−𝜌𝑘4003 + 𝜌𝑙4003 ) Derivace tohoto vztahu podle souřadnic určovaného bodu 4003: −𝑎𝑘𝑙 4003 = −
2.6.3
𝑋40030 − 𝑋 𝑘 𝑋40030 − 𝑋 𝑙 + 𝜌𝑘40030 𝜌𝑙40030
−𝑏𝑘𝑙 4003 = −
𝑌40030 − 𝑌 𝑘 𝑌40030 − 𝑌 𝑙 + 𝑘 𝜌40030 𝜌𝑙40030
−𝑐𝑘𝑙 4003 = −
𝑍40030 − 𝑍 𝑘 𝑍40030 − 𝑍 𝑙 + 𝜌𝑘40030 𝜌𝑙40030
Matice plánu −𝑎14003 −𝑏14003 −𝑐14003 ⎜ ⎟ .. .. .. A=⎝ ⎠, . . . 𝑛 𝑛 𝑛 −𝑎4003 −𝑏4003 −𝑐4003 ⎛
kde 𝑛 je počet 2. diferencí ve všech epochách.
⎞
2.6.4
Opravy přibližných hodnot neznámých (vlastní vyrovnání)
Opravy přibližných hodnot souřadnic přijímače byly vypočteny vyrovnáním MNČ. (︁
dx = − A𝑇 PDD A 2.6.5
)︁−1
· A𝑇 PDD l
Vyrovnané hodnoty souřadnic bodu
Vyrovnané hodnoty neznámých se vypočtou přičtením oprav k přibližným hodnotám. ⎛
⎞
𝑋40030 ⎟ ⎜ x = ⎝ 𝑌40030 ⎠ + dx 𝑍40030 2.6.6
Vektor vyrovnaných měření
S využitím vyrovnaných hodnot neznámých byl vypočten vektor vyrovnaných měření lvyr : 𝑙𝑖𝑣𝑦𝑟 = (𝜌𝑘4004 − 𝜌𝑙4004 ) + (−𝜌𝑘4003 + 𝜌𝑙4003 ), kde 𝜌𝑆𝑅 =
√︁
(𝑋𝑅 − 𝑋 𝑆 )2 + (𝑌𝑅 − 𝑌 𝑆 )2 + (𝑍𝑅 − 𝑍 𝑆 )2 ,
kde 𝑅 jsou stanoviska 4004 a 4003 a 𝑆 je PRN družic. Při výpočtu vyrovnaných měření je nutné zachovat stejné pořadí družic a jako referenční (k) zvolit stejnou družici jako při výpočtu vektoru l. 2.6.7
I. a II. výpočet oprav
Opravy byly vypočteny dvakrát. Jednou v rámci vyrovnání a podruhé jako rozdíl 2. diferencí vypočtených z vyrovnaných neznámých a 2. diferencí měřených pseudovzdáleností. vI = A · dh + l vII = lvyr − lR , kde lR je vektor obsahující 2. diference měřených pseudovzdáleností. 2.6.8
Porovnání I. a II. výpočtu oprav
I. a II. výpočet oprav byly vzájemně porovnány. Jelikož rozdíl řádově 10−6 byl považován za dostatečně malý, nebyla provedena další iterace. vI − vII
2.7
Charakteristiky přesnosti
Na závěr byly vypočteny charakteristiky přesnosti. 2.7.1
Aposteriorní střední jednotková chyba √︃
𝑚0 =
vT v , 𝑛−𝑚
kde n ... počet 2. diferencí (velikost vektoru l) , m=3 ... počet neznámých. K výpočtu byl použit I. výpočet oprav.
2.7.2
Kovarianční matice neznámých (︁
Qdh = AT PDD A 2.7.3
Střední chyby neznámých mdh = 𝑚0 ·
2.8
)︁−1
√︁
𝑑𝑖𝑎𝑔(Qdh )
Číselné výsledky pro bod 4003
Přibližné souřadnice určovaného bodu: 𝑋0 = 3970651.65𝑚 𝑌0 = 1018573.21𝑚 𝑍0 = 4870428.35𝑚 Střední aposteriorní jednotková chyba: 𝑚0 = 0.65 Vyrovnané hodnoty souřadnic přijímače a jejich střední chyby: 𝑋 = 3970658.98𝑚 𝑚𝑋 = 0.41𝑚 𝑌 = 1018588.51𝑚 𝑚𝑌 = 0.27𝑚 𝑍 = 4870437.26𝑚 𝑚𝑍 = 0.46𝑚
3
Shrnutí výsledků
Poloha bodu 4004: 𝑋4004 = 3970663.1𝑚 𝑌4004 = 1018585.9𝑚 𝑍4004 = 4870434.1𝑚
𝑚𝑋 = 2.6𝑚 𝑚𝑌 = 1.7𝑚 𝑚𝑍 = 2.8𝑚
Poloha bodu 4003: 𝑋4003 = 3970658.98𝑚 𝑌4003 = 1018588.51𝑚 𝑍4003 = 4870437.26𝑚
𝑚𝑋 = 0.41𝑚 𝑚𝑌 = 0.27𝑚 𝑚𝑍 = 0.46𝑚
Závěr: V první části úlohy byly z dat obsažených v observačním a navigačním RINEXu určeny souřadnice bodu 4004. K výpočtu byly použity zdrojové kódy pro předchozí dvě úlohy, které byly upraveny tak, aby je bylo možné použít jako funkce. Pro výpočet bylo použito kódové měření a to P-kód na nosné vlně L2. Výsledné souřadnice bodu 4004 jsou určeny s přesností 2-3m. V druhé části úlohy byly určeny souřadnice bodu 4003 a to pomocí dvojitých diferencí. Jako referenční bod byl zvolen bod 4004, jehož souřadnice vypočtené v první části úlohy byly zafixovány. Použitím druhých diferencí byl odstraněn vliv ionosférické refrakce a snížen vliv troposférické refrakce (viz výše). Dále byl eliminován vliv broadcast efemeridů a chyby hodin družice, jejichž střední chyba může být až 1m. Použitím dvojitých diferencí rovněž není potřeba zavádět chybu hodin přijímače jako neznámou do výpočtu, neboť je vyloučena. Při vyrovnání je dosaženo o řád lepší přesnosti výsledných souřadnic. Je však nutné uvážit, že pomocí dvojitých diferencí je určována pouze relativní poloha vůči referenčnímu bodu. Proto i když bylo dosaženo relativní přesnosti v řádu decimetrů, absolutní poloha bodu 4003 je určena se stejnou přesností jako poloha referenčního bodu 4004, tedy 2-3m. To lze demonstrovat porovnáním posunů (𝑑𝑥) obou bodů: ∆𝑋4003 = 7.297𝑚 ∆𝑋4004 = 7.089𝑚 ∆𝑌4003 = 15.308𝑚 ∆𝑌4004 = 15.558𝑚 ∆𝑍4003 = 8.910𝑚 ∆𝑍4004 = 8.764𝑚 Tyto posuny jsou v řádech metrů stejné. Tedy do výsledných souřadnic bodu 4003 se plnou měrou promítly souřadnice bodu 4004. Výpočty byly provedeny v programu Octave. Zdrojový kód k výpočtům není přílohou technické zprávy (v případě potřeby bude zaslán).
V Kralupech nad Vltavou 22.11.2010
Jan Dolista (
[email protected])