+ All Categories
Home > Documents > 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení...

4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení...

Date post: 01-Mar-2021
Category:
Upload: others
View: 3 times
Download: 0 times
Share this document with a friend
15
4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené výsledky za rok 2007. Projekt je rozdělen jak metodicky, tak i časově na tři hlavní fáze: teoretickou, implementační a ověřovací. Vsoučasnosti se řešení nachází ve fázi teoretické, která je plánována dle projektu do 30. 5. 2008, ovšem z důvodů, které nebylo možné řešitelským týmem ovlivnit a které jsou uvedeny podrobněji dále, se navrhuje její prodloužení do 31. 12. 2008. Dále popíšeme podrobněji řešení jednotlivých aktivit vprůběhu roku 2007. Teoretické řešení problematiky geometrických transformací di- gitálních obrazů (Aktivita A07-01) Na základě předchozích výsledků (aktivita A06-01, A06-02) byl navržen výpočetní postup pro stanovení globální charakteristiky přesnosti základních typů geometrických transformací (podobnostní, afinní) a jejich nelineárního zobecnění metodou kolokace. Globální charakteris- tika přesnosti se týká jisté předem dané zájmové oblasti. Formulace problému Vzhledem k tomu, že přesnost transformace bude využita pro vektorové GISy, byl stanoven předpoklad, že zájmová oblast bude zadána jako uzavřený polygon. Za globální charakteristiku přesnosti byla z tradičních důvodů zvolena střední souřadnicová chyba. Střední souřadnicová chyba m XY polohy bodu je obecně definována pomocí jeho kovarianční matice C. m XY := Tr(C) 2 (1) Výpočet hodnoty kovarianční matice C (tj. jejích prvků c ij ) jednoho konkrétního bodu s ohledem na nelineární změny měřítka byl předmětem aktivity A06-02. K výpočtu hodnoty kovarianční matice C je třeba znát příslušné rozdělení pravděpodobnosti onoho konkrétního bodu a jeho střední hodnotu ˆ x. C := R 2 (x - ˆ x) · (x - ˆ x) T f ( x ) dx (2) Symbol f představuje hustotu pravděpodobnosti konkrétního bodu jakožto dvojrozměrného náhodného vektoru.
Transcript
Page 1: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

4.1.1 Popis řešení projektu InGeoCalc za rok 2007

Úvod

V tomto popisu je obsažen postup a dosažené výsledky za rok 2007. Projekt je rozdělen jakmetodicky, tak i časově na tři hlavní fáze: teoretickou, implementační a ověřovací. Vsoučasnostise řešení nachází ve fázi teoretické, která je plánována dle projektu do 30. 5. 2008, ovšemz důvodů, které nebylo možné řešitelským týmem ovlivnit a které jsou uvedeny podrobnějidále, se navrhuje její prodloužení do 31. 12. 2008. Dále popíšeme podrobněji řešení jednotlivýchaktivit vprůběhu roku 2007.

Teoretické řešení problematiky geometrických transformací di-gitálních obrazů (Aktivita A07-01)

Na základě předchozích výsledků (aktivita A06-01, A06-02) byl navržen výpočetní postuppro stanovení globální charakteristiky přesnosti základních typů geometrických transformací(podobnostní, afinní) a jejich nelineárního zobecnění metodou kolokace. Globální charakteris-tika přesnosti se týká jisté předem dané zájmové oblasti.

Formulace problému

Vzhledem k tomu, že přesnost transformace bude využita pro vektorové GISy, byl stanovenpředpoklad, že zájmová oblast bude zadána jako uzavřený polygon. Za globální charakteristikupřesnosti byla z tradičních důvodů zvolena střední souřadnicová chyba. Střední souřadnicováchyba mXY polohy bodu je obecně definována pomocí jeho kovarianční matice C.

mXY :=

√Tr(C)2

(1)

Výpočet hodnoty kovarianční matice C (tj. jejích prvků cij) jednoho konkrétního bodus ohledem na nelineární změny měřítka byl předmětem aktivity A06-02. K výpočtu hodnotykovarianční matice C je třeba znát příslušné rozdělení pravděpodobnosti onoho konkrétníhobodu a jeho střední hodnotu x̂.

C :=∫∫R2(x− x̂) · (x− x̂)T f(x ) dx (2)

Symbol f představuje hustotu pravděpodobnosti konkrétního bodu jakožto dvojrozměrnéhonáhodného vektoru.

Page 2: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

Globální střední souřadnicová chyba má vypovídat o polohové přesnosti nějakého podrob-ného bodu, tj. libovolného bodu uvnitř zájmové oblasti. Proto kovarianční matice C musíbýt pro celou zájmovou oblast stejná. Podle definice (2) však kovarianční matice C závisí napoloze jednoho konkrétního podrobného bodu reprezentované jeho střední polohou x̂ a jemupříslušném rozdělení pravděpodobnosti s hustotou f , které může být v různých místech zájmovéoblasti různé. Proto je třeba místo hustoty pravděpodobnosti f(x) uvažovat podmíněnou hus-totu pravděpodobnosti f(x | t ), která platí pro jeden konkrétní bod t zájmové oblasti. Provýpočet globální střední souřadnicové chyba je tedy nutné nejprve stanovit určité globální roz-dělení pravděpodobnosti pro celou zájmovou oblast. Protože body uvnitř zájmové oblasti jsouspojitě rozloženy, bylo nutno provést plošnou integraci rozdělení pravděpodobnosti přes celouzájmovou oblast.

gT(ε) =∫∫T

f( t+ ε | t )ϕT(t) dt (3)

T . . . daná zájmová oblast

gT . . . globální hustota pravděpodobnosti v oblasti T

ε . . . náhodná chyba bodu v místě t ∈ T,

f( . | t ) . . . podmíněná hustota pravděpodobnosti konkrétního bodu v místě t ∈ T,

ϕT(t) . . . hustota pravděpodobnosti místa t v oblasti T

Hustota pravděpodobnosti ϕT(t) vyjadřuje míru významnosti místa t v rámci oblasti T.Pokud jsou všechny části oblasti T stejně významné, je hustota ϕT na oblasti T rovnoměrná,tj. ϕT(t) = 1

|T| pro všechna t ∈ T, ϕ(t) = 0 pro t 6∈ T . Symbol |T| vyjadřuje plochu oblasti T.Globální rozdělení pravděpodobnosti s hustotou gT zahrnuje jak nepřesnost transformace

v dané zájmové oblasti, tak nejistotu v poloze podrobného bodu modelovanou rozdělením prav-děpodobnosti ϕT. Vzhledem k obecné volbě zájmové oblasti T nemusí být hustota gT v (3) nor-mální. Přesto je možné podle definice (2) vypočítat jí odpovídající globální kovarianční maticiCT.

CT =∫∫R2

ε · εT gT( ε ) dε (4)

Střední souřadnicová chyba mXY vypočtená podle (1 pro kovarianční matici CT je hledanouglobální charakteristikou přesnosti zájmové oblasti T.Velikost globální střední souřadnicové chyby závisí na poloze a velikosti zájmové oblasti.

Vypovídací schopnost globální střední souřadnicové chyby je dána rozložením přesnosti bodův zájmové oblasti. Čím větší jsou rozdíly v přesnosti bodů, tím menší je vypovídací schopnost.Nejlepší vypovídací schopnost má globální střední souřadnicová chyba tehdy, když přesnostbodů je všude v zájmové oblasti přibližně stejná.

Page 3: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

Řešení problému

Těžiště výzkumné práce spočívalo ve výpočtu dvojných integrálů (3), (4). Vzhledem k nepřílišvysokému počtu dimenzí obou integrálů by bylo možné moužít některou z univerzálních me-tod numerické integrace, avšak díky některým příznivým okolnostem to nebylo nutné. Jednouz těchto okolností je skutečnost zjištěná v průběhu aktivity A0602, že rozdělení pravděpodob-nosti f( . | t ) je přibližně normální. Pokud je navíc rozložení přesnosti v zájmové oblasti dosta-tečně rovnoměrné, bude přibližně normální i globální rozdělení pravděpodobnosti gT. V takovémpřípadě lze i rozdělení pravděpodobnosti gT aproximovat normálním rozdělením, případně roz-vojem v dvojrozměrnou Edgeworthovu řadu (viz [2], kapitola 3).Výpočet integrálu (3) probíhá za pomoci polynomické aproximace změny hustoty pravděpo-

dobnosti f( . |t ) při měnícím se t ∈ T. Polynom 3. stupně aproximuje prvky kovarianční maticerozdělení f( . | t ). Za těchto zjednodušujících předpokladů byly odvozeny přibližné vzorce provýpočet plošného integrálu nad polygonální oblastí, které umožnily řešit integrál (3) analy-ticky. Integrál (4) pak bylo možno vyčíslit pomocí dvojrozměrné Edgeworthovy řady a tu pakintegrovat člen po členu. Součástí přibližného řešení je i odhad chyby aproximace.

Závěr

Globální střední souřadnicová chyba v polygonální oblasti T se vypočte postupně podle vzorců(3), (4), (1). K výpočtu je použita speciální metoda numerické integrace založená na aproximaciintegrandu hustotou normálního rozdělení pravděpodobnosti.

Teoretické řešení problematiky klasifikace digitálních obrazů(Aktivita A07-02)

Celková koncepce projektu si vynutila netradiční přístup ke klasifikaci digitálních obrazů.Existuje nepřeberné množství rozmanitých metod klasifikace. Zvolená jednotná metodika, zalo-žená na bayesovském přístupu, přirozeně vedla k volbě bayesovské klasifikace. Obvyklý postupbayesovské klasifikace však bylo nutno poněkud modifikovat. Byl vynechán jeden z finálníchkroků tohoto postupu, a sice optimalizace průběhu hranic jednotlivých klasifikovaných oblastí.Výsledkem optimalizace jsou totiž jednoznačně stanovené hranice oblastí a tato jednoznačnostneumožňuje odhadnout nepřesnost průběhu hranic.

Formulace problému

Hranice klasifikovaných oblastí jsou hledány ve tvaru uzavřených polygonů, jejichž vrcholynemají jednoznačně stanovenou polohu. Poloha těchto vrcholů je dána rozdělením pravděpo-dobnosti.

Page 4: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

Obrázek 1: Polygonální oblast s nepřesnou hranicí

Úkolem modifikované bayesovské klasifikace je tedy nalézt rozdělení pravděpodobnosti vr-cholů všech polygonů aproximujících příslušné oblasti v digitálním obraze. Tato úloha je prin-cipiálně velmi složitá, a proto její součástí musí být i volba vhodných zjednodušujících předpo-kladů, které umožní danou úlohu efektivně řešit.

Řešení problému

Hlavním metodickým nástrojem znázorňování oblasti s nepřesnou hranicí se stal pojem fuzzymnožina. Protože mezivýsledkem bayesovské klasifikace je tzv. pravděpodobnostní mapa, (ma-tice pixelů ohodnocených pravděpodobností, že příslušný pixel náleží k určité třídě). Tato prav-děpodobnostní mapa byla použita jako tzv. funkce příslušnosti k fuzzy množině. To znamená,že stupněm příslušnosti k fuzzy množině je pravděpodobnost. Takováto pravděpodobnostníinterpretace fuzzy množin bývá sice zastánci teorie fuzzy množin většinou odmítána, avšakv průběhu našeho výzkumu se toto skloubení teorie fuzzy množin a teorie pravděpodobnostiukázalo být velmi přínosným a užitečným.Nejprve byl navržen rámcový algoritmus, jak na základě nepřesné hranice oblasti sestrojit

fuzzy množinu, tj. funkci příslušnosti jednotlivých pixelů k dané třídě odpovídající uvažovanéoblasti. Tento algoritmus je založen na velmi složitém vzorci obsahujícím mnohonásobný ur-čitý integrál jisté funkce, jejíž analytický tvar není jednoznačně definován. Hodnoty tohotointegrandu však musí odpovídat rozdělení pravděpodobnosti polohy vrcholů polygonu aproxi-mujícího uvažovanou oblast. Aby bylo možno řešit zmíněný mnohonásobný integrál analyticky,byly stanoveny následující zjednodušující předpoklady:

Page 5: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

1. Rozdělení pravděpodobnosti vrcholů polygonů jsou isotropní neboli kruhová, tzn. majístejný rozptyl ve všech směrech.

2. Hustotou pravděpodobnosti vrcholů je polynomiální funkce odpovídající Pearsonovu roz-dělení typu II.

3. Topologie rozložení oblastí v obraze (tj. sousednost jednotlivých oblastí nezávisle na jejichtvaru) je známa předem, např. na základě předběžné, přibližné klasifikace některou zestandardních metod.

Bylo prokázáno, že za těchto předpokladů je možno zmíněný mnohonásobný integrál analytickyvyřešit a sestrojit tak fuzzy množiny odpovídající jednotlivým oblastem.Dalším krokem bylo řešení inverzní úlohy k předchozímu kroku, a sice nalezení rozdělení

pravděpodobnosti vrcholů polygonů na základě znalosti odpovídajících funkcí příslušnosti. Tatoinverzní úloha je opět řešena pomocí bayesovského přístupu (viz např. [1]). Odhadovanými pa-rametry jsou v tomto případě koeficienty Pearsonova rozdělení. Data představují pravděpodob-nostní mapy získané v prvním kroku bayesovské klasifikace.Podstatu této inverzní úlohy lze nejlépe znázornit na příkladě jednorozměrné oblasti, což

ukazuje obrázek 2.

0 2 4 6 8 10 120.0

0.2

0.4

0.6

0.8

1.0

Obrázek 2: Odhad hranice oblasti

Pravděpodobnostní mapu hledané oblasti představuje žlutý histogram v pozadí. Červenoulomenou čarou je zakreslena schodovitá funkce příslušnosti k fuzzy množině, která má mode-lovat hledanou oblast. Této funkci příslušnosti odpovídají dva krajní body nepřesné hraniceznázorněné modrou a černou zvonovitou křivkou. Změnou tvaru červené schodovité funkce pří-slušnosti se mění poloha a strmost obou krajních zvonovitých křivek. Přitom se mění také

Page 6: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

pravděpodobnost, že při dané fuzzy množině naměříme data reprezentovaná žlutým histogra-mem. Bayesova věta nám umožní vypočíst inverzní pravděpodobnosti, tj. pravděpodobnostirozmanitých tvarů červené schodovité funkce příslušnosti za předpokladu pevně daných dat.Protože mezi tvarem červené schodovité funkce příslušnosti a krajními zvonovitými křivkamiexistuje jednoznačná závislost, lze takto určit rozdělení pravděpodobnosti polohy zvonovitýchkřivek i jejich rozptylu. Tím je odhadnuta hranice hledané oblasti i její nepřesnost. V reálnémdvojrozměrném případě je princip řešení inverzní úlohy stejný, pouze související výpočty jsoumnohem složitější.Matematický zápis právě popsaného řešení byl publikován v [3].Aby bylo možno považovat hodnoty pravděpodobností z pravděpodobnostních map za ná-

hodné veličiny, byly stanoveny další zjednodušující předpoklady:

4. Hodnoty pravděpodobnostních map jsou považovány za relativní četnosti.

5. Rozdělení pravděpodobnosti absolutních četností je multinomické rozdělení.

6. Statistická závislost sousedních pixelů je zanedbatelná.

Za těchto předpokladů je možno odhadnout parametry rozdělení pravděpodobnosti jednot-livých vrcholů polygonů, avšak zatím pouze s pomocí numerické maximalizace celkové hustotypravděpodobnosti odhadovaných parametrů.

Závěr

Navržená modifikace bayesovské klasifikace zahrnuje dvojí použití Bayesovy věty. Nejprve sepoužije diskrétní verze Bayesovy věty (tzv. Bayesův vzorec) k výpočtu pravděpodobnostníchmap hledaných oblastí. Podruhé se Bayesova věta uplatní při odhadu polohy vrcholů polygo-nální hranice tak, aby se odpovídající fuzzy množina co nejlépe přimykala k pravděpodobnostnímapě vypočtené v předchozím kroku (jako červená schodovitá funkce ke žlutému histogramuna obrázku 2).Zmíněné zjednodušující předpoklady, které usnadňují výpočet za cenu snížení obecnosti,

nejsou definitivní, v průběhu dalšího výzkumu budou hledány způsoby, jak je oslabit. Je protopotřebné prodloužit dobu výzkumu v rámci dílčího cíle V001 (matematická formulace postupuvýpočtů).

Reference

[1] Karl-Rudolf Koch. Bayesian Inference with Geodetic Applications, volume 31 of LectureNotes in Earth Sciences. Springer–Verlag, 1990.

[2] K. V. Mardia. Families of Bivariate Distributions, volume 27 of Griffin’s Statistical Mono-graphs & Courses. Griffin, London, 1970.

[3] L. Soukup. Overall uncertainty of georeferencing and classification. In A. Stein, editor,International Symposium on Spatial Data Quality. ISSDQ ’07, Enschede, June 2007. ITCEnschede.

Page 7: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

Teoretické řešení problematiky transformace výšek a určení vybraných parametrů tíhového pole Země (Aktivita A07-03) Tato aktivita navazovala na aktivitu A06-04. Výsledkem je seznam parametrů tíhového pole Země, které bude možné pomocí znalostního systému určovat z dostupných geodat uložených v databázi znalostního systému.

Vybrané parametry jsou uvedeny v tab. 1.

parametr výpočet z lok.

modelu výpočet z

glob. modelu W tíhový potenciál NE ANO U potenciál normálního pole NE ANO T poruchový potenciál NE ANO g velikost tíhového zrychlení ANO NE η příčná složka tížnicové odchylky ANO NE ξ meridiánová složka tížnicové odchylky ANO NE Δg tíhová anomálie ANO ANO gr radiální derivace tíhového zrychlení NE ANO AT terénní korekce ANO NE AB gravitační zrychlení Bouguerovy sfér. slupky ANO NE ΔgB Bouguerova anomálie ANO NE ζ výšková anomálie NE ANO ζ převýšení kvazigeoidu ANO NE N převýšení geoidu ANO NE

Tabulka 1. Vybrané parametry tíhového pole Země

V příloze 1 jsou uvedeny konkrétní postupy, jak určovat jednotlivé parametry z

lokálních a globálních modelů. U každého z parametrů bylo nutno rozhodnout, zda se bude počítat z globálního

modelu nebo z lokálního modelu. Globální modely jsou k dispozici ve formě Stokesových koeficientů, současný model EGM96 (s koeficienty do stupně a řádu 360, tj. prostorové rozlišení 30’x30’) měl být v průběhu roku 2007 nahrazen mnohem přesnějším modelem EGM07 (s koeficienty do stupně a řádu 2160, tj. s prostorovým rozlišením 5’x5’). Vydání nového modelu EGM však bylo odloženo na rok 2008, nyní je již nový model s označením EGM08 hotov a bude představen na EGU General Assembly ve Vídni v dubnu 2008. V červnu 2008 na sympoziu GGEO v Chanii na Krétě bude jedna sekce věnována právě testování tohoto nového modelu. Přestože se naše pracoviště podílelo na testování nového globálního modelu, publikace jakýchkoliv výsledků před jeho oficiálním představením na konferenci ve Vídni je zakázána. Nicméně pro účely této zprávy lze potvrdit to, co bylo možné očekávat: přesnost nového modelu výrazně překračuje přesnost modelu EGM96. Model EGM96 tedy bude v tomto projektu v příštím roce plně nahrazen novým modelem EGM08. Na území ČR dosahuje nový globální model takové přesnosti v hodnotách tíhové anomálie, že dokonce překračuje přesnost některých starších lokálních modelů. Z tohoto důvodu bude pro výpočet většiny parametrů ve znalostním systému používán právě tento model. Pouze pro parametry, kde je nezbytná vysoká přesnost, nebo kde nelze globální modely použít, se bude vycházet z lokálních modelů. Jedná se např. o hodnoty terénní korekce AT , které vyjadřují gravitační vliv topografických hmot a počítají se z digitálních

Page 8: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

modelů terénu a modelů laterálních variací hustot topografických hmot. Dále se jedná o velikost tíhového zrychlení g (a další z g přímo odvozované veličiny jako např. tíhové anomálie) v případě, že budou požadovány hodnoty na topografii, protože je přesnější je interpolovat z databáze středních hodnot g, která vychází z terestrických gravimetrických měření. Konečně se jedná o parametry: převýšení geoidu N, převýšení kvazigeoidu ζ a tížnicové odchylky η, ξ, kde je klíčová přesnost, a proto je stále vhodnější vycházet z předpočítaných lokálních modelů. Zejména hodnota ζ, která bude využívána pro nejdůležitější funkci systému, kterou je transformace výšek mezi elipsoidálními a normálními (Moloděnského) výškami. Shrnutí Pro úlohu transformace výšek a výpočet parametrů tíhového pole Země je pro účely projektu popsán matematický aparát výpočtu jednotlivých veličin z globálních nebo lokálních modelů. Zároveň jsou v databázi uloženy přesné modely geoidu, kvazigeoidu, tížnicových odchylek, tíhového zrychlení, výšek, laterálních variací hustoty topografických hmot, terénní korekce, globální model tíhového pole Země a další parametry uvedené v příloze 1. Cíle aktivity byly splněny. V současnosti probíhá poměrně rychlý rozvoj jak globálních, tak lokálních modelů, a proto budou jednotlivé modely v dalších fázích řešení projektu dále průběžně aktualizovány. Aktualizace modelů bude nutná i po uvedení aplikace do provozu. Odstranění systematické složky chyby modelu kvazigeoidu a odhad chyby transformace výšek (Aktivita A07-04) Při využívání lokálního modelu gravimetrického kvazigeoidu pro transformaci elipsoidálních výšek (měřených např. pomocí GPS) je nutno odstranit systematický posun. Tento systematický posun je způsoben různými posuny referenčních hladin národních výškových systémů proti referenční hladině tíhového systému, která nemůže být přesně určena, protože nelze měřit přímo hodnotu potenciálu gravitačního pole, ale pouze velikost jeho gradientu (velikost tíhového zrychlení).

V plánu na tento rok jsme měli pokusit se o využití Bayesovského přístupu k určení tohoto systematického posunu. Pokusy o aplikaci Bayesovského odhadování neznámých parametrů však narazily na neznalost rozdělení pravděpodobnosti odhadované veličiny. Vzhledem k tomu, že o tomto posunutí nemáme žádnou apriorní znalost, lze sice teoreticky Bayesovský odhad použít, ovšem ve skutečnosti by nevedl ke zpřesnění odhadu, ale naopak by zavedl do odhadu neznámého parametru podmínky, jež by neměly fyzikální opodstatnění. To souvisí i s neznalostí statistických charakteristik vstupních tíhových dat – databáze tíhových dat, kterou máme k dispozici, obsahuje pouze velikost tíhového zrychlení v rastru, který byl získán interpolací z bodových měřených hodnot zemské tíže. Interpolace tíhových dat je ovšem velice komplexní úloha, protože hodnoty tíže silně závisí na výšce a nelze je proto interpolovat přímo, ale nejprve se musí numericky odstranit závislost měřených hodnot na výšce, následně provést interpolaci a pak zpětně závislost na výšce pro jednotlivé buňky rastru zrekonstruovat. K tomu se použijí střední hodnoty výšek digitálního modelu terénu v rastru, který svým rozlišením odpovídá rozlišení generovaného rastru tíže.

Odhad statistických vlastností takto interpolovaných středních hodnot tíže je velice obtížný i za předpokladu, že bychom znali přesnou metodu tvorby rastru tíhových dat (použitý způsob redukce tíhových měření, způsob interpolace, atd.). My bohužel přesnou metodiku tvorby neznáme, a přesnost tíhových dat tedy nedovedeme apriori odhadnout.

Page 9: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

Odhad neznámého parametru jsme tedy provedli pomocí střední hodnoty rozdílu mezi výškami kvazigeoidu určenými gravimetricky a mezi výškami kvazigeoidu určenými na bodech výběrové údržby, na nichž jsou známy jak elipsoidální výšky v systému WGS84, tak normální Moloděnského výšky změřené přesnou nivelací. Postup odstranění systematické chyby lokálního modelu kvazigeoidu Varianta 1 Pro testování bylo využito 1024 bodů výběrové údržby, na nichž je známa elipsoidální výška nad elipsoidem WGS84 (měřeno GPS) a nivelovaná výška v systému Balt po vyrovnání (Bpv). Převýšení kvazigeoidu z nich určíme jednoduše podle vztahu

nivelWGSGPS Hh −=ζ Tyto hodnoty budeme porovnávat s hodnotami ζgrav z přesného lokálního modelu kvazigeoidu spočítaného P. Novákem v r. 2007. Abychom mohli provést srovnání, z lokálního modelu kvazigeoidu ve formě rastru vyinterpolujeme hodnoty ζgrav v bodech výběrové údržby. Protože model kvazigeoidu má poměrně vysoké rozlišení, je rastr poměrně spojitý a nezáleží příliš na volbě interpolační metody (tj. na počtu okolních buněk rastru, z nichž interpolujeme výšku ζgrav)

interpolace

lineární nejbližší soused kubická

střední hodnota -0.3891 -0.3888 -0.3890 směrodatná odchylka 0.0566 0.0574 0.0566 minimum -0.8911 -0.8850 -0.8911 maximum 0.5570 0.5730 0.5571

Tab. 1: Statistika rozdílů ζGPS - ζgrav na 1024 bodech výběrové údržby [m]

-1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.60

50

100

150

200

250

300

Obr. 3: Histogram rozdílů ζGPS - ζgrav na 1024 bodech výběrové údržby

Střední hodnota rozdílů ζGPS - ζgrav vyjadřuje systematický posun výšek ζGPS a ζgrav. Ze statistik a histogramu je zřejmé, že soubor rozdílů obsahuje odlehlé hodnoty, které bude nutno odstranit.

Page 10: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

Odstranění odlehlých pozorování Pro detekci a odstranění odlehlých hodnot byly použity 2 testy: iterační Grubbsův test s hladinou významnosti 5 % a 3 σ test (test, který odstraňuje všechny hodnoty ve vzdálenosti větší než trojnásobek směrodatné odchylky od střední hodnoty). Statistiky po odstranění těchto hodnot jsou uvedeny v tabulce 2 (pro test 3 σ byla využita hodnota σ po odstranění odlehlých pozorování pomocí Grubbsova testu). Výskyt bodových odlehlých měření lze přisoudit buďto možným chybám určení GPS výšky nebo změně výšky bodu vybrané údržby (Určení výšky bodu pomocí nivelace a GPS se od sebe mohlo časově podstatně lišit).

Grubbsův test test 3σ

střední hodnota -0.3888 -0.3897směrodatná odchylka 0.0431 0.0415minimum -0.5123 -0.5123maximum -0.2197 -0.2599odstraněných hodnot 4 10

Tab. 2: Statistika rozdílů ζGPS - ζgrav po odstranění odlehlých hodnot [m]

Obr. 4: Znázornění detekovaných odlehlých hodnot, vlevo Grubbsovým testem,

vpravo testem 3σ Po odstranění odlehlých pozorování lze odhadnout velikost systematické chybu modelu ζ0 = -0,39 m. Směrodatná odchylka charakterizuje přesnost transformace GPS výšek do národního systému Bpv. Střední chyba transformace výšek by byla přibližně 4 cm. Po odstranění odlehlých pozorování pomocí testu 3σ byla testována normalita souboru rozdílů. Na základě testu Kolmogorov-Smirnov nelze na hladině významnosti 5 % vyvrátit hypotézu o normalitě souboru rozdílů, zatímco hypotéza o normalitě souboru rozdílů s odlehlými pozorováními byla na téže hladině významnosti zamítnuta. Všech 10 detekovaných odlehlých pozorování jsme tedy odstranili a dále pracujeme pouze s 1014 body výběrové údržby. Pokud vykreslíme všechny hodnoty, které nepadnou do intervalu ohraničeného ±2,5násobkem směrodatné odchylky od střední hodnoty (viz obr. 3), ukáže se, že největší chyby jsou koncentrovány v okolí hranic ČR. To bylo možné očekávat ze dvou důvodů:

1. Na hranicích ČR se nachází hory, kde je průběh kvazigeoidu komplikovanější. 2. Při výpočtu lokálního modelu jsme v příhraničních oblastech využívali i data

z okolních států, která nemáme v takové přesnosti jako na území ČR.

Page 11: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

Obr. 5: Body vybrané údržby,na kterých | ζGPS - ζgrav | nepadne do intervalu ohraničeného ±2,5násobkem směrodatné odchylky od střední hodnoty

Uvážíme-li přesnost elipsoidálních výšek, byť měřených GPS statickou, systematické chyby výšek určených nivelací a chyby v tíhových datech (zejm. chyby způsobené jejich interpolací do pravidelného rastru), lze považovat střední chybu transformace kolem 4 cm za uspokojivou. Varianta 2 Důvod pro odstranění systematického posunu mezi ζgrav a ζGPS má zřejmé fyzikální opodstatnění (viz výše). Pro transformaci výšek měřených pomocí GPS do Bpv je však vhodné řešení gravimetrického geoidu upravit tak, aby lépe odpovídalo naměřeným hodnotám ζGPS. Obr. 4 zobrazuje residuální rozdíly mezi hodnotami (ζgrav+ ζ0) a ζGPS (po odstranění odlehlých pozorování). Z obrázku je zřejmé, že residuální rozdíly mají na částech území ČR systematický charakter. Tato systematická posunutí mohou být způsobena např. lokálními deformacemi výškového systému Bpv a budou pro účely transformace výšek odstraněny.

Obr. 6: Residuální rozdíly mezi hodnotami ζGPS a (ζgrav+ ζ0) (pro větší názornost byly bodové

hodnoty interpolovány do rastru) Místo určování konstantního posunu mezi ζgrav a ζGPS (viz předchozí odstavce) přičítáme k výškám ζgrav korekční kubický polynomiální člen závislý na souřadnicích ϕ, λ

Page 12: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

00

3

0

3

0,,),( λλλϕϕϕλϕλϕδζ −=Δ−=ΔΔΔ=∑∑

= =n m

nmmna

Tímto způsobem byl určován i průběh kvazigeoidu CR2000 (viz [1]). Souřadnice ϕ0, λ0 jsou zvoleny přibližně uprostřed území ČR (ϕ0=50˚, λ0=15˚) a koeficienty amn jsou určeny vyrovnáním rozdílů ζgrav a ζGPS pomocí MNČ na bodech vybrané údržby.

Postup vede ke snížení rozdílů mezi ζgrav a ζGPS (předchozí případ byl vlastně variantou tohoto způsobu, kde ovšem všechny koeficienty kromě a00 byly rovny 0) a rovněž směrodatná odchylka rozdílů mezi (ζgrav+δζ(ϕ, λ)) a ζGPS bude menší v předchozím případě.

Koeficient a00 (opět představuje posunutí výšek nezávisle na poloze) vychází a00=0.3887, tedy opět velice podobně jako ve variantě 1, všechny ostatní koeficienty vychází v řádu setin nebo tisícin. Statistika rozdílů mezi (ζgrav+δζ) a ζGPS je uvedena v tabulce 3.

řád interpol. polynomu

0 (konstantní posunutí) ζGPS - (ζgrav+ζ0)

3 (kubický polynom) ζGPS - (ζgrav+δζ)

střední hodnota 0 0 směrodatná odchylka 0.0415 0.0214 minimum -0.1226 -0.0989 maximum 0.1298 0.0985

Tab. 3: Statistika rozdílů ζgrav a ζGPS po zavedení opravných členů [m]

Obr. 7: Residuální rozdíly mezi hodnotami ζGPS a (ζgrav+ δζ) (pro větší názornost byly bodové hodnoty interpolovány do rastru)

Závěr Na základě relativně přesného lokálního modelu gravimetrického kvazigeoidu vypočteného P. Novákem byl odvozen systematický posun mezi gravimetrickým geoidem a tíhovými anomáliemi vypočtenými z rozdílu nivelace a GPS měření. Pro transformaci výšek byla dále odvozena varianta modelu, který je upraven tak, aby lépe odpovídal známým tíhovým anomáliím ζGPS. Tato varianta modelu umožňuje transformovat elipsoidální výšky určené pomocí GPS do systému Bpv se střední chybou 2 cm. Cíle aktivity byly částečně splněné,

Page 13: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

nebyla splněna část týkající se aplikace Bayesovského odhadování parametrů. Důvody spočívají v neznalosti statistických vlastností tíhových dat. Vzhledem k provedenému rozsáhlému testování modelu kvazigeoidu pomocí nezávislých dat lze považovat střední chybu transformace výšek i odhad systematického posunu za spolehlivě určené. Poznámka k přesnosti tíhových dat Při testech dostupných tíhových dat (databáze VÚGTK), které prováděl prof. P. Novák s prof. R. Kleesem v r. 2007 se však ukázalo, že tíhová data uložená v modelu VÚGTK jsou méně přesná (odhad střední chyby vyšel v řádu mGal), než se očekává u dat získaných z pozemních tíhových měření. Protože tíhová data byla měřena s vysokou přesností, chyba by mohla být způsobena jejich interpolací. Z těchto tíhových dat byl počítán tento, ale i některé další modely lokálního gravimetrického kvazigeoidu na území ČR. Pokud by se podařilo zpřesnit databázi pozemních tíhových měření, lze očekávat ještě přesnější model gravimetrického kvazigeoidu. Vývoj lokálního modelu kvazigeoidu bude tedy dále pokračovat.

Řešený projekt bude výšky transformovat z předem připraveného rastru, který bude samozřejmě možné v průběhu řešení projektu i po jeho spuštění průběžně aktualizovat tak, aby zahrnoval vždy aktuální modely kvazigeoidu pro ČR. Reference [1] J. Kostelecký, J. Kostelecký Jr., I. Pešek, J. Šimek, O. Švábenský, J. Weigel, A. Zeman: Quasigeoids for the territory of the Czech Republic. Studia Geophysica et Geodaetica, č. 48, 3/2004, s. 503-518, ISSN 0039-3169.

Teoretické řešení problematiky analýzy posunů a deformací. (Aktivita A07-05) Cílem této aktivity bylo převést zjištěné posuny bodů na souvislé pole deformací. Doposud používaný způsob interpolace ve čtvercové síti, jak je uvedeno například v [1] a na obrázcích 8 a 9, není vhodný pro odhad přesnosti posunu v libovolně zvoleném (obecném) bodě pole deformací, neboť výpočet posunu v obecném bodě by musel být zprostředkovaný výpočtem posunů v bodech čtvercové sítě. Navíc by se interpolace musela provádět nadvakrát: jednou pro body čtvercové a z nich pak podruhé pro obecné body sítě.

V rámci jednotné metodiky projektu je vhodnější použít metodu kolokace, jejíž parametry jsou bayesovsky odhadovány. Při tom se uplatní některé teoretické poznatky získané v aktivitách A06-01, A06-02. Posuny v poli deformací lze interpretovat jako transformaci souřadnic. Praktické řešení takovéto transformace je však mnohem složitější než v aktivitách A06-01, A06-02, neboť v tomto případě je transformace silně nelineární. Proto také velmi záleží na správném určení statistické závislosti sousedních bodů. Nelze tudíž použít obvyklý postup stanovení kovarianční matice založený na jedné neměnné kovarianční funkci. Tento obvyklý postup totiž předpokládá homogenní a isotropní průběh vzájemné závislosti blízkých bodů v poli deformací. Tento předpoklad nemusí být při deformaci splněn. Míra statistické závislosti v poli deformací je ovlivňována nejen polohovou blízkostí příslušných bodů, ale též velikostí posunu v daném místě. Aby byla zohledněna tato

Page 14: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

skutečnost, byl navržen původní postup stanovení kovarianční matice, při němž je kovarianční funkce proměnná v závislosti na poloze obecného bodu i na směru spojnice ovlivňujících se bodů.

Obr. 8: Příklad pole posunů interpolovaných do čtvercové sítě

Obr. 9: Příklad pole deformací určených pro čtvercovou síť

Page 15: 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 · 2011. 7. 28. · 4.1.1 Popis řešení projektu InGeoCalc za rok 2007 Úvod V tomto popisu je obsažen postup a dosažené

Takto modifikovaná metoda kolokace ponechává větší míru uplatnění naměřeným datům (posunům) a potlačuje vliv celkového trendu. Celkový trend posunu nebyl odhadován. Buďto je předem znám s dostatečnou přesností (např. z předchozích měření posunů) anebo je dán fyzikálními podmínkami experimentu. Odhadovanými parametry jsou tedy pouze koeficienty kovarianční funkce. Bylo testováno několik typů kovarianční funkce. Vzhledem k nelinearitě kovarianční funkce bylo nutné při testování používat numerickou integraci. Nakonec byla vybrán polynomiální typ kovarianční funkce, jejíž koeficienty se nejsnáze odhadují. Reference [1] Talich, M.: Geometrical Analysis of Deformation Measurement using Continuum Mechanics by Web Application. In: Strategic Integration of Surveying Services, The XXX FIG General Assembly and Working Week 2007, Hong Kong SAR, China, 13-17 May 2007, ISBN 978-87-90907-59-4, http://www.fig.net/pub/fig2007/papers/ts_1f/ts01f_03_talich_1375.pdf

Závěr Projekt se v r. 2007 stále nacházel v teoretické fázi svého řešení. Vzhledem k nevydání nového přesnějšího modelu tíhového pole Země, který měl být vydán již v roce 2007 pod názvem EGM 2007, ale bude vydán až v dubnu 2008 jako EGM08, se navrhuje posunout termín dokončení dílčího cíle V001 "matematicky formulovat postup výpočtů" o půl roku na 31. 12. 2008.

V roce 2007 byl navržen výpočetní postup pro stanovení globální charakteristiky přesnosti základních typů geometrických transformací (podobnostní, afinní) a jejich nelineárního zobecnění metodou kolokace. Za globální charakteristiku přesnosti byla z tradičních důvodů zvolena střední souřadnicová chyba, která vypovídá o polohové přesnosti nějakého podrobného bodu, tj. libovolného bodu uvnitř zájmové oblasti. V oblasti teoretického řešení problematiky klasifikace digitálních obrazů se využívá bayesovské klasifikace a byl navržen originální postup stanovení hranic klasifikovaných oblastí, který respektuje jejich neurčitost.

Teoretické řešení problematiky transformace výšek a určení vybraných parametrů tíhového pole Země pokračovalo odvozením vzorců pro výpočty některých dalších parametrů tíhového pole Země podle zadání projektu, které nebyly počítány v r. 2006.

Probíhaly práce na odstranění systematické složky chyby modelu kvazigeoidu a odhadu chyby transformace výšek. Bylo dosaženo přesnosti při transformaci elipsoidálních výšek určených pomocí GPS do systému Bpv se střední chybou 2 cm.

Při teoretickém řešení problematiky analýzy posunů a deformací bylo dosaženo vyjádření průběhu pole deformací spojitou funkcí, takže je možno přesně určit parametry tenzoru deformace v libovolném bodě zájmové oblasti.

Ve všech oblastech řešení projektu jsou odvozovány potřebné výpočetní vzorce tak, aby bylo možné navrhnout výpočetní algoritmy a postupy vedoucí k softwarovému zabezpečení řešení a realizaci vlastního znalostního systému.


Recommended