ČESKÁ ZEMĚDĚLSKÁ UNIVERZITA V PRAZE
Fakulta životního prostředí
_________________________________________________________________________
Diplomová práce
Hydrologicky korektní model terénu povodí Modrava 1 na základě trojúhelníkové sítě Petr Herout _________________________________________________________________________
Vedoucí práce: Ing. Jiří Pavlásek, Ph.D.
Obor: Environmentální modelování Duben 2013
i
ii
Prohlášení Prohlašuji, že jsem diplomovou práci na téma „Hydrologicky korektní model terénu povodí Modrava 1 na základě trojúhelníkové sítě“ vypracoval samostatně za použití uvedené literatury.
V Praze dne 21. dubna 2013 ……………………… Petr Herout
iii
Diploma Thesis
Terrain model of Modrava 1 catchment based on triangular network
Presented MSc. diploma thesis is focused on Digital Terrain Models (DTM) based on the Triangular Irregular Network (TIN) and on the hydrological and terrain analysis of these models. In the theoretical part topics of the digital terrain modeling and of the triangulation used for creating TIN are reviewed. In the second (practical) part a model of the experimental basin of Modrava 1 (Sumava mountin, Southern Bohemia) and a set of used algorithms are introduced. The set of algorithms could be applied on the creating of similar models of small catchments. In the model there are flow paths of the surface water, shape, area, slope and orientation relations of the catchment analyzed. The study is focused on the evaluation of the model quality loss in dependence on the terrain generalization rate in the end.
Key words: Triangular Irregular Network – TIN – Digital Terrain Model – Modrava 1 Catchment – Triangulation
iv
Poděkování Tímto bych rád poděkoval Ing. Jiřímu Pavláskovi, PhD., vedoucímu mé diplomové práce, za vstřícný přístup, za jeho ochotu a trpělivost. Dále chci poděkovat celému týmu, který se mnou podílel na terénním měření, zejména pak Petru Baštovi a Kubovi Medkovi, kteří práci v terénu obětovali mnoho hodin. Panu Ing. Danielu Zahradníkovi, PhD. patří dík za inspiraci a probuzení mého zájmu o programování. A v neposlední řadě chci poděkovat rodině a přátelům, kteří mi byli oporou nejen po celou dobu tvorby této diplomové práce. Děkuji vám P.H.
v
1
ÚVOD
-2-
1.1
CÍLE DIPLOMOVÉ PRÁCE
-2-
1.2
STRUKTURA PRÁCE
-3-
TEORETICKÁ ČÁST
-4-
2.1
ZÁKLADNÍ POJMY
-4-
2.2
DIGITÁLNÍ MODELY TERÉNU
-6-
2
2.2.1 2.3
DATOVÉ SÍTĚ MODELŮ TRIANGULACE
-7- 11 -
2.3.1
FORMULACE PROBLÉMU A ZÁKLADNÍ VLASTNOSTI TRIANGULACÍ
- 11 -
2.3.2
ROZDĚLENÍ TRIANGULAČNÍCH ALGORITMŮ
- 12 -
3
METODIKA
- 19 -
3.1
ZÁJMOVÁ OBLAST
- 20 -
3.2
SBĚR DAT
- 22 -
3.3
PŘÍPRAVA DAT
- 23 -
3.4
STRUKTURA PROGRAMU
- 24 -
3.5
VYBRANÉ ALGORITMY A VÝPOČTY
- 26 -
3.5.1
ODSTRANĚNÍ BEZODTOKÝCH OBLASTÍ
- 26 -
3.5.2
VYTVOŘENÍ A ANALÝZA SEZNAMU TROJÚHELNÍKU
- 27 -
3.5.3
ANALÝZA HRAN
- 29 -
3.5.4
NALEZENÍ TRASY ODTOKU
- 31 -
3.5.5
NALEZENÍ ROZVODNICE
- 32 -
3.5.6
VÝPOČET PLOCHY POVODÍ
- 33 -
3.5.7
VERIFIKACE MODELU
- 34 -
PRÁCE S PROGRAMEM
- 35 -
3.6.1
PŘÍPRAVA PROSTŘEDÍ
- 35 -
3.6.2
TVORBA MODELU TERÉNU
- 36 -
3.6.3
ANALÝZA MODELU
- 37 -
3.6.4
INTERAKTIVNÍ FUNKCE
- 38 -
3.6
vi
3.6.5
ULOŽENÍ A NAČTENÍ MODELU
- 39 -
3.6.6
FUNKCE PRO OVLÁDÁNÍ GRAFIKY
- 40 -
3.6.7
PŘÍKLAD POUŽITÍ
- 41 -
4
VÝSLEDKY A VÝSTUPY DIGITÁLNÍHO MODELU TERÉNU
4.1
VÝSTUPY ANALÝZY TERÉNU
- 42 - 43 -
4.1.1
BODOVÉ POLE
- 43 -
4.1.2
TROJÚHELNÍKOVÁ SÍŤ
- 44 -
4.1.3
ODSTRANĚNÍ BEZODTOKÝCH OBLASTÍ
- 45 -
4.1.4
ORIENTACE TROJÚHELNÍKŮ
- 46 -
4.1.5
SKLONY TROJÚHELNÍKŮ
- 48 -
4.1.6
ANALÝZA HRAN
- 50 -
4.2
VÝSTUPY ANALÝZY POVODÍ
- 50 -
4.2.1
ROZVODNICE
- 51 -
4.2.2
TVAROVÉ CHARAKTERISTIKY POVODÍ
- 52 -
4.3
VÝSLEDKY VERIFIKACE A POROVNÁNÍ MODELŮ
- 56 -
4.3.1
POROVNÁNÍ MODELŮ VYTVOŘENÝCH NA BÁZI TROJÚHELNÍKOVÉ SÍTĚ
- 59 -
4.3.2
SROVNÁVACÍ MODEL NA BÁZI ČTVERCOVÉ SÍTĚ
- 60 -
5
DISKUZE
- 61 -
5.1
MĚŘICKÉ METODY
- 61 -
5.2
DIGITÁLNÍ MODEL TERÉNU
- 62 -
5.3
ANALÝZA POVODÍ
- 63 -
6
ZÁVĚR
- 65 -
7
SEZNAM INFORMAČNÍCH ZDROJŮ
- 66 -
8
SEZNAM PŘÍLOH
- 70 -
vii
1 Úvod Výběr tématu mé diplomové práce ovlivnily tři hlavní faktory. Prvním z nich byla láska k přírodě, která mě vedla k rozhodnutí, že při tvorbě diplomové práce nesmím strávit veškerý čas u počítače a v knihovnách, nýbrž alespoň část musím strávit nějakou prací v terénu. S tím úzce souvisí druhý faktor mého rozhodnutí. Již v posledním ročníku bakalářského studia jsem poprvé zavítal na jedno z šumavských experimentálních povodí, která provozuje Katedra vodního hospodářství a environmentálního modelování Fakulty životního prostředí České zemědělské univerzity. Tehdy jsem zde pomáhal při infiltračních pokusech a krátce také při geodetickém měření a již v té době se v mé hlavě zrodila myšlenka, že by se moje diplomová práce mohla týkat některého z těchto povodí. Třetím faktorem byla záliba v hlavolamech a rébusech, díky které se jedním z mých nejoblíbenějších předmětů v prvním ročníku studia Environmentálního modelování stal předmět Počítačové modelování. Během jednoho roku jsem si díky tomuto předmětu a hlavně díky velmi pozitivnímu přístupu vyučujícího, pana Ing. Daniela Zahradníka, PhD., okruh mých zájmů rozšířil o programování. Později, když jsem se dověděl, že povodí Modrava 2 je jediné ze čtyř šumavských povodí KVHEM, které je zmapováno geodetickým měřením s dostatečnou podrobností, a je také jediné, pro které byl vytvořen podrobný digitální model terénu, rozhodl jsem se, že právě tvorba digitálního modelu terénu některého ze zbývajících experimentálních povodí bude tématem mé diplomové práce. O důvodu výběru povodí Modrava 1 bude pojednáno později v kapitole Zájmová oblast. Vzhledem k tomu, že nerad chodím po vyšlapaných cestách, jsem se rozhodl, že místo abych pro tvorbu modelu užil komerční software, raději vytvořit vlastní počítačový program.
1.1 Cíle diplomové práce Hlavním cílem práce je geodetické zaměření povodí Modrava 1 a tvorba digitálního modelu terénu povodí na bázi nepravidelné trojúhelníkové sítě pomocí
-2-
vlastního souboru algoritmů. Podružnými cíli je vyšetření směrů odtoku povrchových vod, nalezení orografické rozvodnice a zjištění některých terénních a hydrologických ukazatelů povodí; to vše na základě vytvořeného digitálního modelu. Vytvořený program by měl být univerzální, použitelný pro libovolnou sadu výškových dat. Vedlejším, osobním cílem je rozvoj mých programovacích znalostí a schopností práce s programem R.
1.2 Struktura práce Text práce je rozdělen na několik základních celků. Kapitola 2 Teoretická část se v první podkapitole zabývá shrnutím problematiky tvorby digitálních modelů terénu s největším zaměřením na modely na bázi nepravidelné trojúhelníkové sítě. Další podkapitola se zaměřuje na triangulační algoritmy a jejich výhody a nevýhody. V obou případech je uvedena pouze stručná rešerše s odkazy na další doporučenou literaturu.
Třetí kapitola Metodika podrobně popisuje postup prací na tvorbě digitálního modelu experimentálního povodí Modrava 1 od sběru dat, přes přípravu dat pro další použití, až po popis jednotlivých použitých algoritmů. Dále je v této kapitole uveden podrobný návod pro práci s předkládaným programem.
Ve čtvrté kapitole jsou shrnuty výsledky a výstupy digitálního modelu terénu. Vzhledem k tomu, že výstupy mají převážně formu objemných datových tabulek, které jsou v textové formě pouze nepřehlednou změtí čísel, jsou výstupy znázorňovány povětšinou pouze formou grafů, diagramů a obrázků. Právě tyto obsáhlé obrazové výstupy jsou také důvodem značného rozsahu práce.
-3-
2 Teoretická část Obsahem této kapitoly je shrnutí základních pojmů používaných v tomto textu a dále úvod do problematiky a stručná rešerše k tématům digitální modely terénu a triangulační algoritmy.
2.1 Základní pojmy Povodí Hrádek a Kuřík (2002) popisují povodí jako základní hydrologickou oblast, ve které je zkoumán hydrologický proces, a ve které je zjišťován vzájemný vztah bilančních prvků, a definují ho jako území z hydrologického hlediska uzavřené, do kterého nepřitéká žádná voda po povrchu ani pod povrchem. Je ohraničené rozvodnicí. Veškeré srážky, které na povodí spadnou, se povrchovým i podpovrchovým odtokem dostávají do tzv. uzavírajícího (nebo také uzávěrového) profilu, kde povodí opouštějí. Uzávěrový profil je nejnižším bodem celého povodí (Hrádek – Kuřík 2002).
Rozvodnice Rozvodnice je myšlená čára, která značí hranici mezi povodími. Rozlišujeme rozvodnici orografickou, která určuje povodí povrchových vod, a rozvodnici hydrologeologickou, která ohraničuje povodí podpovrchových vod. Orografická rozvodnice počíná v uzávěrovém profilu a probíhá terénem po hřbetnicích, hřebenech, vrcholech a sedlech. Průběh hydrogeologické rozvodnice je závislý na uložení nepropustných vrstev a na geologické stavbě území (Hrádek – Kuřík 2002). Hydrogelogické a orografické povodí náležící ke stejnému uzávěrovému profilu se mohou v menší či větší míře lišit. Pokud je dále v této práci řeč o povodí či o rozvodnici, je vždy míněno povodí povrchových vod a rozvodnice orografická.
-4-
Plocha povodí Plocha je nejzákladnější geometrickou charakteristikou povodí. Jedná se o plochu průmětu povodí do vodorovné roviny. Nejčastěji se udává v jednotkách [km2] a stanovuje se planimetrováním z mapy (Moore – Grayson – Ladson 1991, Hrádek – Kuřík 2002). Vzhledem k malé ploše zájmového území je však v této práci udávána v [m2]. Způsob vypočtení plochy povodí v této práci je popsán v kapitole 3 Metodika.
Délka povodí Délkou
povodí
se
nazývá
vzdálenost
od
nejvyššího
bodu
povodí
k uzávěrovému profilu (Moore – Grayson – Ladson 1991).
Sklon svahu, orientace svahu Sklon svahu v bodě je definován jako vertikální odchylka tečné plochy terénu v daném bodě od vodorovné roviny. Lze ji uvádět v procentech nebo ve stupních. Vhledem k nejednotnosti tvaru povodí se sklonové poměry vyjadřují charakteristikou střední sklon svahů v povodí, která se často určuje jednoduchým vztahem:
I SV =
H max − H min F
⋅ 100 [%],
(1)
kde Hmax a Hmin je maximální a minimální nadmořská výška v povodí [m] a F je plocha povodí [m2] (Hrádek – Kuřík 2002). Další možností je využití Herbstova vzorce, který počítá s průměrnou délkou vrstevnic:
I SV =
∆h ⋅ ∑ l si F
⋅ 100 [%],
(2)
kde ∆h je výškový interval mezi vrstevnicemi, lsi průměrná délka vrstevnic v itém intervalu a F je plocha povodí (Hrádek – Kuřík 2002).
-5-
Orientací svahu je míněna expozice svahu ke světovým stranám. V literatuře se lze setkat i s označením aspekt (angl. aspect) (Moore – Grayson – Ladson 1991, Bašta 2008). V bodě ji lze nejlépe definovat jako pravostrannou odchylku průmětu spádnice v daném bodě do vodorovné roviny od severu (neboli azimut).
2.2 Digitální modely terénu Názvosloví týkající se digitálních modelů reprezentujících terén není v literatuře zcela jednotné. Kromě výrazu digitální model terénu (Digital terrain model, DMT) se můžeme setkat například s výrazy digitální model reliéfu, digitální výškový (elevační)
model
(Digital
elevation
model,
DEM),
digitální
model
povrchu.
Terminologický slovník zeměměřictví a katastru nemovitostí Výzkumného ústavu geodetického, topografického a kartografického (VÚGTK 2013) tyto názvy vykládá doslova takto:
Digitální model terénu / reliéfu (DMT,DMR) (Digital terrain model, DTM): 1: digitální reprezentace zemského povrchu v paměti počítače, složená z dat a interpolačního algoritmu, který umožňuje mj. odvozovat výšky mezilehlých bodů 2: v USA má zpravidla formu nepravidelně rozmístěných výškových bodů, které optimálně vystihují terénní tvary včetně hran a výškových extrémů
Digitální výškový model (Digital elevation model, DEM): 1: digitální model reliéfu pracující výhradně s nadmořskými výškami bodů 2: datová sada výškových hodnot, které jsou algoritmicky přiřazeny k 2rozměrným souřadnicím 3: v USA soubor nadmořských výšek ve vrcholech mříže vytvořené v pravidelných intervalech souřadnic x a y národního referenčního souřadnicového systému
Digitální model povrchu (Digital surface model, DSM):
-6-
zvláštní případ digitálního modelu reliéfu konstruovaného zpravidla s využitím automatických prostředků (např. obrazové korelace ve fotogrammetrii) tak, že zobrazuje povrch terénu a vrchní plochy všech objektů na něm (střechy, koruny stromů apod.).
Podle tohoto lze říci, že anglicky psaná literatura rozlišuje rozdíl mezi výrazy Digital terrain model (DTM) a Digital elevation model (DEM) jako rozdíl v datové struktuře modelu (viz níže). V praxi se však setkáváme s vnímáním digitálního modelu terénu jako uspořádaného číselného pole představujícího prostorové rozdělení některých charakteristik terénu. Digitální elevační model je pak výraz podřazený předchozímu (DTM) a popisuje případ, kdy je jako charakteristika terénu použita výhradně nadmořská výška terénu (Moore – Grayson – Ladson 1991, Wilson – Gallant 2000). V některé anglické literatuře je pak zkratka DTM užívána pro výraz Digital terrain modeling, což je obecný výraz pro veškeré techniky používané pro tvorbu digitálních elevačních modelů. Výrazy digitální model terénu a digitální elevační model tato literatura vnímá jako synonyma (Hengel – Gruber – Shrestha 2003 in Bašta 2008).
2.2.1 Datové sítě modelů Digitální elevační modely můžeme rozdělit podle struktury datové sítě do tří základních typů (Obrázek 1).
Obrázek 1: Struktury datových sítí: a) pravidelná síť, b) nepravidelná trojúhelníková síť a c) síť na bázi vrstevnic (Moore – Grayson – Ladson 1991)
Nejužívanějšími z nich jsou modely na bázi pravidelné sítě (regular grid). Výšková data jsou zde strukturována do pravidelné sítě diskrétních bodů. Sítě mohou
-7-
nabývat různých tvarů, nejčastěji se jedná o sítě čtvercové, alternativně však lze použít i sítě obdélníkové, dále čtyřhranné neortogonální či trojúhelníkové a hexagonální, přičemž může být v každém případě volena libovolná hustota datové sítě. Pro každý bod pravidelné sítě jsou uloženy hodnoty kartézských souřadnic x, y a z (Bašta 2008). Důvodem velkého rozšíření tohoto typu struktury modelu je především velice snadná implementace výpočetních algoritmů a také relativně malá výpočetní náročnost. Nevýhodami pravidelně strukturovaných sítí pak jsou především špatná schopnost zachycovat náhlé změny v elevaci terénu, ovlivnění výsledků rozlišením sítě (stejně tak i výrazné ovlivnění výpočetní náročnosti), nerealistické výsledky některých algoritmů hydrologických analýz a v neposlední řadě jednotnost rozlišení sítě v celé rozloze modelovaného území (Moore – Grayson – Ladson 1991). V českém jazyce jsou možnosti použití modelů na bázi pravidelné sítě a algoritmy hydrologických analýz na nich používaných velmi dobře shrnuty v pracích Bašty (2008) a Bartáka (2008). Další
velmi
rozšířenou
datovou
strukturou
DEM
je
nepravidelná
trojúhelníková síť (Triangulated irregular network, TIN). Používá ji například Tajchman (1981), Jones et al. (1990) a Yu et al. (1997). Struktura dat TIN je založena na bodech (tzv. uzlové body, anglicky nodes), které jsou nepravidelně rozmístěné na povrchu modelovaného terénu. Všechny tyto uzlové body jsou pomyslnými úsečkami spojeny do sítě, která celý povrch modelovaného terénu dělí v nepravidelné trojúhelníky (v anglicky psané literatuře zvané facets). Aby bylo dosaženo optimálního výsledku je vhodné, aby uzlové body byly umístěny v charakteristických místech terénu jako jsou: vrcholy, sedla, hřebeny, údolí (údolnice), lomové hrany atp. (Peucker – Fowler – Little 1978, Moore – Grayson – Ladson 1991). Metody používané pro tvorbu trojúhelníkové sítě (tzv. triangulace) budou podrobněji rozebrány v kapitole 2.3. Srovnání metod pro extrakci bodů pro TIN modely z modelů na bázi pravidelného gridu provedl Lee (1991). Na rozdíl od pravidelného gridu jsou v TIN modelech kromě x, y a z kartézských souřadnic ukládány ještě odkazy na sousední trojúhelníky. Ve výsledku je tedy třeba uložit tři datové tabulky: v první jsou zapsány kartézské souřadnice uzlových bodů, ve druhé jsou pro každý trojúhelník ukazatele na tři uzlové body a ve třetí jsou pro každý trojúhelník odkazy na tři sousední trojúhelníky (pokud leží trojúhelník na okraji modelu, je jedno místo odkazu prázdné). Existuje i druhá alternativa uložení, ve které je jako základní
-8-
kámen brán uzlový bod namísto trojúhelníku. V takovém případě jsou pak pro každý bod uzlový uloženy kromě kartézských souřadnic odkazy na všechny jeho sousední body (připojené úsečkami) v daném pořadí (např. od nejsevernějšího souseda po směru hodinových ručiček). Druhá varianta je datově méně náročná, avšak je náročnější na implementaci výpočetních algoritmů hydrologických analýz (Peucker – Fowler – Little 1978). Souřadnice z bývá mezi uzlovými body nečastěji lineárně interpolována z vrcholů – uzlových bodů – trojúhelníku, kterému náleží. Lze však použít i specifické metody počítající se sklonem sousedních trojúhelníků, např. Bezierovy pláty Coonsovy plochy (Kilimánek 2006). Na bázi TIN modelů lze provádět nepřeberné množství analýz terénu. Nejběžnějšími analýzami jsou automatické stínování map (hill shading), analýza sklonu svahů (slope mapping), tvorba vrstevnic (contouring), analýza viditelnosti (line of sight) (Peucker – Fowler – Little 1978) a z hydrologických analýz pak například můžeme jmenovat vyhledání rozvodnice (watershed derive), nebo trasování odtoku (derive streem lines) (Vivoni 2004, Kilimánek 2006). Výhodou modelů na bázi TIN oproti modelům na bázi pravidelné sítě je především schopnost popisovat povrch v různých úrovních rozlišení – v místech s větší členitostí terénu může být síť bodů libovolně zhuštěna a naopak; to zapříčiňuje větší efektivitu v ukládání dat. Navíc vzhledem k umístění uzlových bodů do význačných bodů terénu dochází u TIN modelů k lepší reprezentaci terénu než u gridu. Nevýhodou je složitější implementace výpočetních technik (Moore – Grayson – Ladson 1991). Třetí používanou datovou strukturou je síť na bázi vrstevnic (contour-based network). Ta se skládá z vrstevnic digitalizovaných do podoby tzv. Digital line graphs (DLG) a ze spádnic. Spádnice a vrstevnice jsou navzájem ortogonální a terén rozdělují na síť nepravidelných polygonů. Tuto strukturu poprvé představili Onstad a Brakensiek (1968, in Moore – Grayson – Ladson 1991). Její výhodou je především sledování přirozených tras odtoku, díky čemuž je možné zjednodušit složité třídimenzionální rovnice proudění na soustavu vzájemně závislých jednodimenzionálních rovnic. Proto se tato struktura uplatňuje hlavně v některých hydrologických aplikacích (Wilson – Gallant 2002, Moore – O’Loughlin – Burch 1988).
-9-
Hydrologicky korektní model terénu můžeme definovat jako model, který neobsahuje žádné terénní deprese (angl. sink), ať už vzniklé nepřesnou reprezentací terénu nebo takové , které se v terénu opravdu nachází. V modelech na bázi gridu často vznikají deprese vinou chyb při interpolaci dat. Depresí se u nich rozumí buňka, jejíž všechny sousední buňky mají vyšší hodnotu elevace (dno deprese), případně skupina buňek, jejichž okolní buňky splňují stejnou podmínku (uzavřená deprese, bezodtoká oblast) (Barták 2008). Analogicky se dá uvažovat o depresích a bezodtokých oblastech i v modelech na bázi TIN. Místo sousedních buněk stačí posuzovat sousední uzlové body modelu. Pokud model bezodtoké oblasti obsahuje, je potřeba je před prováděním všech hydrologických analýz odstranit.
- 10 -
2.3 Triangulace S pojmem triangulace se můžeme setkat v mnoha různých významech napříč nejrůznějšími vědními obory. V tomto textu je však triangulací míněno rozdělení plochy nebo povrchu na trojúhelníky dle určitých pravidel. Triangulace je používána v nejrůznějších aplikacích například v kartografii (tvorba digitálních modelů terénu), v počítačové grafice (vizualizace prostorových dat), v biometrii (detekce otisků prstů), ke zpracování obrazu (segmentace a rozpoznávání vzorů) a v mnoha dalších aplikacích. V této kapitole poskytnu pouze stručný úvod do problematiky, která je velmi obsáhlá. Podrobnější informace je nutné dohledat v dalších zdrojích, kterých je k dispozici například v internetových knihovnách nepřeberné množství. Valná většina zdrojů se však věnuje výhradně jednomu konkrétnímu algoritmu. V českém jazyce obecně problematiku triangulačních algoritmů shrnuje velmi dobře Kolingerová (1999 a 2003) a Zábranský (2005). Z anglicky psané literatury více druhů algoritmů popisují Hjelle a Dæhlen (2010), ti navíc poskytují široký přehled další vhodné literatury. Stručně ale velmi přehledně popisuje většinu základních algoritmů také Shojaee (et al. 2006).
2.3.1 Formulace problému a základní vlastnosti triangulací Máme dánu množinu bodů P={p1,p2,…,pn} v R2 a hledáme triangulaci T nad množinou P. Triangulace T nad množinou bodů P představuje takové planární rozdělení, které vytvoří soubor trojúhelníků t = {t1,t2,…,tm} tak, aby platila všechna následující pravidla. (1) Libovolné dva trojúhelníky ti, tj ∈ T, ti ≠ tj mají společnou nejvýše jednu stranu. (2) Sjednocení všech trojúhelníků t∈ T tvoří v euklidovském prostoru souvislou množinu. (3) Uvnitř žádného trojúhelníku neleží další bod P (Bayer 2013, Hjelle – Dæhlen 2010). Počet hran a počet trojúhelníků konvexní triangulace je závislý na počtu uzlových bodů. Řídí se vztahy: N t = 2N − k − 2
(3)
a
- 11 -
N h = 3N − k − 3 ,
(4)
kde Nt značí počet vzniklých trojúhelníků, Nh počet vzniklých hran, N počet uzlových bodů a k počet bodů tvořících konvexní obálku (Peucker – Fowler – Little 1978, Bayer 2013).
2.3.2 Rozdělení triangulačních algoritmů Triangulační algoritmy patří mezi nejvíce teoreticky rozpracované postupy, přesto však vznikají stále nové algoritmy, které se co nejlépe snaží splňovat požadavky, které jsou na ně uplatňovány. Nejdůležitějšími požadavky na triangulace jsou především jednoduchost algoritmu, snadná implementace, vysoká výpočetní rychlost i pro značně objemné sady vstupních dat, malá citlivost na singulární případy a optimální tvar trojúhelníkové sítě. Některé z těchto požadavků stojí v kontrastu: jednoduchý algoritmus se snadnou implementací může jen těžko dosahovat vysokých výpočetních rychlostí (Bayer 2013). Triangulace můžeme rozdělovat z hlediska kvality různých kritérií na (1) lokálně optimální, (2) globálně optimální a (3) multikriteriálně optimální. Pro lokálně optimální triangulace platí, že každý čtyřúhelník tvořený dvojicí sousedních trojúhelníků je optimalizován vzhledem k vybranému kritériu. Kritérii pro lokální optimalizaci nejčastěji bývají minimální nebo maximální úhel v trojúhelníku, minimální nebo maximální výška trojúhelníku, poloměr opsané kružnice, poloměr vepsané kružnice, plocha trojúhelníku atd. Ve dvojici trojúhelníků může tedy být například maximalizován minimální úhel. Pro globálně optimální triangulace platí, že zadané kritérium musí být optimalizováno pro všechny trojúhelníky v celé triangulaci. Mezi globální kritéria patří například suma délek všech hran. Někdy je za globální kritérium považováno také kritérium obsahu tzv. povinných hran, tedy hran, které nesplňují lokální kritérium, avšak jsou povinně zahrnuty do triangulace. Multikriteriálně optimální triangulace kombinují více globálních či lokálních kritérií. Nejsou pro ně doposud známé dostatečně efektivní algoritmy (Bayer 2013, Kolingerová 1999). Běžněji jsou však v literatuře triangulace děleny podle jejich geometrické konstrukce. Mezi nejznámější patří (1) hltavá triangulace (greedy triangulation), (2)
- 12 -
Delaunayho
triangulace,
(3)
triangulace
s minimální
váhou
(Minimum
Weight
Triangulation, MWT), (4) triangulace s povinnými hranami (Constrained Triangulation) a (5) datově závislé triangulace (Data Dependent Triangulation, DDT).
2.3.2.1 Hltavá triangulace Hltavá neboli greedy triangulace (GT) je složená z nejkratších možných, vzájemně se neprotínajících hran. Snadno se implementuje, avšak výpočetní složitost bývá velká, nicméně při vhodné úpravě algoritmu lze ji výrazně snížit (Dickerson – Drysdale – McElfresh 1998). U GT není brán zřetel na tvar trojúhelníků, a proto se příliš nehodí pro tvorbu digitálních modelů terénu. Výhodou je, že se do ní snadno vkládají povinné hrany (viz kapitola 2.3.2.4). Pokud žádné dvě hrany nemají stejnou délku, je triangulace jednoznačná (Zábranský 2005). Algoritmus GT je jednoduchý a přímočarý; žádná z hran, které se postupně začleňují do triangulace, není dále měněna. Vypadá následovně:
1. Vytvoř seznam všech potencionálních hran a seřaď jej vzestupně podle jejich délky. 2. Vyber první hranu v seznamu (nejkratší) a pokud se nekříží s žádnou již přijatou hranou, zařaď ji do triangulace. 3. Opakuj krok 2, dokud není vybrán potřebný počet hran.
Obrázek 2: Hltavá triangulace. Na obrázcích je patrné postupné včleňování hran. (Bayer 2013)
- 13 -
2.3.2.2 Delaunayho triangulace Delaunayho triangulace (DT) je v současné době jednou z nejpoužívanějších triangulací v nejrůznějších odvětvích včetně digitálního modelování terénu. Poprvé byla představena Borisem Nikolaevichem Delaunayem (anglická transkripce ruského jména Делоне) již v roce 1934. DT lokálně maximalizuje minimální úhel v trojúhelnících, čímž je dosaženo ideálního, neprotáhlého tvaru trojúhelníků. Pro danou množinu bodů je DT vždy jednoznačná, pokud žádné čtyři body neleží na jedné kružnici. Jedno ze základních pravidel DT je, že pro každý trojúhelník musí platit, že uvnitř kružnice jemu opsané nesmí ležet žádný další bod triangulace. Okraje DT jsou vždy shodné s její konvexní obálkou (Zábranský 2005, Shojaee – Helali – Alesheikh 2006).
Obrázek 3: Delaunayho triangulace se znázorněnými opsanými kružnicemi. (Bayer 2013)
Pro tvorbu Delaunyaho triangulací lze použít více různých metod. (1) Prohazování stran, (2) metoda rozděl a panuj, (3) inkrementální vkládání, (4) inkrementální konstrukce a další. Metoda prohazování stran (Edge flip) se zakládá na libovolné triangulaci, kterou lokálně vylepšuje až do konečné podoby DT. Jsou posuzovány čtyřúhelníky složené ze dvou sousedních trojúhelníků. V případě, že trojúhelníky nesplňují zadané kritérium, je
- 14 -
ve čtyřúhelníku použita opačná úhlopříčka. K vylepšování může být užito pravidlo o opsané kružnici případně lokální kritérium minimálního úhlu (Kolingerová 1999).
Obrázek 4: Princip prohazování stran (edge flip) u nevyhovujících trojúhelníků. (Bayer 2013)
Metoda Rozděl a panuj (Devide-and-Conquer) množinu vstupních bodů rozděluje na dvě přibližně stejně veliké podmnožiny, ty dále dělí a takto rekurzivně pokračuje, až zbude v každé podmnožině dostatečně malé množství bodů. Nad každou z těchto základních podmnožin vytvoří DT složenou pouze z několika trojúhelníků. Poté množiny bodů a triangulace opět spojuje a v místech spojů opravuje lokální nedostatky prohazováním hran na základě lokálního kritéria minimálního úhlu trojúhelníků. Výhodou tohoto algoritmu je nízká výpočetní náročnost (Hjelle – Dæhlen 2010).
Metoda inkrementálního vkládání v prvním kroku vytváří obalový trojúhelník, který obsáhne všechny zadané body. Následně do triangulace postupně vkládá náhodně zvolené body ze zadané množiny. Po vložení bodu nejprve vyhledá, ve kterém trojúhelníku současné triangulace leží. Poté přicházejí v úvahu dvě varianty: (1) Pokud leží uvnitř trojúhelníku, vzniknou tři nové trojúhelníky spojením aktuálního vkládaného bodu se všemi vrcholy trojúhelníku, ve kterém leží. (2) Pokud leží přesně v hraně triangulace, dva trojúhelníky, které tuto hranu sdílí, jsou rozděleny hranami vedoucími od vkládaného bodu do protějšího vrcholu trojúhelníku. Poté jsou všechny dotčené trojúhelníky
- 15 -
legalizovány případným prohozením stran pomocí pravidla o prázdné opsané kružnici. Pokud dojde u některého trojúhelníku k prohození, musí být provedena legalizace i u jeho sousedů. Takto pokračuje vkládání jednotlivých bodů, dokud nejsou všechny obsaženy v triangulaci. Na závěr je už jen odstraněn obalový trojúhelník a DT je hotová (Hjelle – Dæhlen 2010, Zábranský 2005).
Obrázek 5: Inkrementální vkládání. Na horní dvojici obrázků je znázorněn případ, kdy vkládaný bod je uvnitř trojúhelníku, na dolní dvojici je případ, kdy se vkládaný bod nachází na hraně. (Bayer 2013)
Na počátku metody inkrementální konstrukce je vybrán libovolný bod ze zadané množiny. Často je jako pivot vybírán bod s extrémními souřadnicemi x nebo y. K němu je nalezen nejbližší bod, se kterým vytvoří směrově orientovanou úsečku a zároveň první hranu triangulace. Následně je na pravé straně od orientované úsečky hledán bod, který splňuje pravidlo o prázdné opsané kružnici, pokud není nalezen, změní se směr orientované úsečky a totéž se opakuje na opačné straně. Když je nalezen odpovídající bod, je spojen s aktuálně zkoumanou úsečkou do trojúhelníku, který je zahrnut do triangulace. Nově vzniklé hrany dostávají souhlasnou orientaci s první orientovanou úsečkou a jsou zařazeny do fronty. Následně je z fronty vždy odebrána první úsečka, na které se opakuje proces hledání dalšího vhodného bodu. Nově vzniklé hrany jsou opět řazeny na konec fronty. To se opakuje dokud, není fronta prázdná (Kolingerová 2003).
- 16 -
Obrázek 6: Metoda inkrementální konstrukce. Případ pivota s extrémní souřadnicí. (Bayer 2013)
Existuje nepřeberné množství metod tvorby DT, výše zmíněné jsou však nejběžnější. Některé další příklady uvádí například Leifer (2006).
2.3.2.3 Triangulace s minimální váhou Triangulace s minimální váhou (Minimal Weight Triangulation, MWT) je taková triangulace zadané množiny bodů, jejíž konvexní obal je shodný s konvexním obalem zadané množiny a součet délek všech hran triangulace je nejmenší možný. Lloyd (1977) dokázal, že ani Delaunayho triangulace ani greedy triangulace zadané podmínky MWT nesplňují. Tvorba MWT je velice náročná a dosud patří k ne zcela vyřešeným problémům výpočetní geometrie. Většinou se zakládá na výchozí greedy triangulaci, která je sice lokálně optimální vzhledem ke kritériu minimálních délek stran, avšak globálně optimální není (Zábranský 2005).
2.3.2.4 Triangulace s povinnými hranami Triangulace s povinnými hranami (Constrained triangulation) není vlastně samostatnou kategorií. Jedná se o modifikace jiných typů – nejčastěji Delaunayho
- 17 -
triangulace s povinnými hranami (Delaunay Constrained Triangulation, DCT) a greedy triangulace s povinnými hranami (Constrained Greedy Triangulation, CGT). Jsou to triangulace, které respektují některé předem zvolené hrany, pro které musí být splněna podmínka, že se nesmí vzájemně křížit. V některých případech povinné hrany umožňují věrněji reprezentovat modelovaný povrch. Vstupem těchto triangulací je kromě množiny bodů ještě množina takových hran. Tyto hrany modelují významné prvky terénní kostry, jako jsou hřbetnice a údolnice, a geomorfologické tvary – například terénní zlomy, které mohou jinak být běžnou triangulací zanedbány a deformovány. Tvorba CGT je snadná. Povinné hrany jsou do triangulace zahrnuty již před začátkem algoritmu a nově zahrnované hrany se jim podřizují. Algoritmy DCT jsou na implementaci o poznání složitější. Lze například nejdříve vytvořit Delaunayho triangulaci a poté v místech povinných hran provádět lokální prohazování. Jejich použití popisují například Gudmundsson (et al. 2005) a Hjelle a Dæhlen (2010).
2.3.2.5 Datově závislé triangulace Všechny dosud zmíněné triangulace byly striktně planární, vstupovaly do nich body zadané pouze kartézskými souřadnicemi x a y. Delaunayho triangulace sice vytváří vhodně tvarované trojúhelníky pro modelování terénu, ovšem právě pouze v rovině. Když se triangulovaným bodům následně přiřadí hodnota elevace, dochází k deformaci trojúhelníků a tím ke změně jejich vnitřních úhlů. V případě, že se v modelovaném terénu nenacházejí trojúhelníky s velkým sklonem, je tato deformace zanedbatelná. V opačném případě, dochází k výraznějším deformacím a tvar trojúhelníků přestává být ideální. Datově závislé triangulace (Data Dependent Triangulation, DDT) obsahují informace o výšce a podřizují jim optimalizační algoritmy. Vstupem bývá většinou již hotová triangulace (vzhledem k dobrému tvaru sítě to bývá Delaunayho triangulace) doplněná o zmíněná výšková data. Poté je tato triangulace optimalizována pomocí specializovaných kritérií, které popisují například Hjelle a Dæhlen (2010) nebo Zábranský (2005). Výhodou DDT je, že nemusí být vkládány povinné hrany, aby model více odpovídal realitě, výrazné body terénní kostry jsou algoritmem automaticky detekovány.
- 18 -
3 Metodika Pro zhotovení digitálního modelu terénu byl vytvořen vlastní soubor skriptů a příkazů spustitelný ve statistickém programovacím prostředí R. Tento soubor je v textu práce pro zjednodušení nazýván programem, ačkoliv se o žádný samostatně spustitelný program nejedná. V první své podkapitole se kapitola Metodika zabývá informacemi o zájmovém území, v dalších pak už postupy a algoritmy, kterými bylo dosaženo výsledného modelu terénu. Následuje podrobný návod na užívání předkládaného programu.
Model terénu byl nejprve vytvořen z kompletní sady změřených dat, následně byly jejím zredukováním na 75, 50 a 25 procent vytvořeny tři další zdrojové soubory dat. Z každého z nich byl opět vytvořen jeden model. Základní model je dále nazýván Model M1, ostatní v závislosti na procentuálním podílu použitých bodů Model M1_75, Model M1_50 a Model M1_25. Všechny čtyři modely byly verifikovány a porovnány s úmyslem zjistit, do jaké míry by bylo možné terén zájmového povodí dále generalizovat, aniž by došlo k výraznější ztrátě přesnosti modelu.
Pro srovnání dvou typů modelů (TIN a grid) byl vytvořen také jednoduchý model na bázi pravidelné čtvercové sítě. Ten byl tvořen v komerčním softwaru ArcGIS 10 pomocí několika specializovaných funkcí. Aby bylo možné tento i postup zopakovat, uvádím v příloze 1 názvy použitých funkcí. Samotné funkce ovšem pro celou práci nestačí, je proto potřeba i základních znalostí a dovedností v programu ArcGIS.
- 19 -
3.1 Zájmová oblast Zájmová oblast byla volena mezi třemi eventualitami. Povodí Modrava 2 bylo už zaměřeno dříve (Bašta 2008), zbývala tedy dvě modravská povodí a povodí Pastouška. Povodí Modrava 1 bylo vybráno z důvodu, že v mrtvém lese, který toto povodí pokrývá, začíná mladý smrkový podrost dosahovat většího vzrůstu a byla tedy poslední příležitost, zaměřit povrch pomocí totální stanice, než to bude podrostem znemožněno úplně. Na zbývajících povodích se nachází vzrostlý les, a proto je možnost měření touto metodou omezená. Povodí Modrava 1 bylo – spolu s dalšími dvěma modravskými povodími – založeno Katedrou vodního hospodářství a Katedrou biotechnických úprav krajiny FLE ČZU roku 1998 v rámci výzkumných aktivit grantového projektu VaV 620/6/97 „Obnova biodiverzity a stability lesních ekosystémů v pásmu přirozeného rozšíření smrku na území NP Šumava“ (KVHEM 2013). Účelem bylo sledovat odlišnosti komponent srážkoodtokového procesu na územích s různým pokryvem. Lokalita Modrava 1 byla silně postižena kůrovcovou kalamitou, byla zde vyhlášena bezzásahová zóna NP Šumava, zůstal zde stát odumřelý smrkový porost, který nyní postupně podléhá větrným kalamitám. V povodí označovaném Modrava 2, byla povolena těžba kalamitního dřeva a zůstala zde paseka, která byla posléze zalesněna smrkem s příměsí klenu a jeřábu. Lokalita Modrava 3 nebyla zasažena kůrovcovou kalamitou vůbec, nachází se zde přibližně stopadesátiletý smrkový porost s příměsí buku (Pavlásek – Máca – Ředinová 2006). Zájmové povodí Modrava 1 se nachází v pramenné oblasti Roklanského potoka (hydrologické pořadí 1-08-01-006) (Pavlásek – Máca – Ředinová 2006), na samé hranici České republiky s německým Bavorskem, asi tři kilometry východně od vrcholu Velkého Roklanu (1453 m n.m.) a dva kilometry jižně od Medvědí hory (1224 m n.m.). Lokalita leží na severním úbočí bezejmenného vrchu ležícího mezi Velkým Roklanem a Blatným vrchem (Vojenský kartografický ústav 1991).
- 20 -
Obrázek 7: Lokalizace zájmového území na mapě. Zdroj: mapy.cz
Část suchých kmenů odumřelého lesa již vlivem větru popadala; probíhá zde postupná obnova smrkového lesa, některé mladé stromy již dosahují i větší než pětimetrové výšky. V horní části povodí se nachází také jednotlivě přimíšené jeřáby. Povrch terénu je pokryt travním porostem a popadanými suchými kmeny a větvemi. Hloubka půdního profilu je 0,4 – 0,6 m, půdním typem je podzol s nadložními horizonty s vyšší mírou zrašelinění (Pavlásek – Máca – Ředinová 2006). Dosud uváděné charakteristiky povodí byly odvozeny manuálně ze základní mapy 1:10 000, pro výpočet středního sklonu svahů byl použit vzorec dle Herbsta (Pavlásek 2013).
Tabulka 1: Charakteristiky povodí zjištěné manuálním odečtením z mapy. KVHEM 2013
Plocha povodí
0,10 km2
Min. nadmořská výška
1216 m n.m.
Max. nadmořská výška
1270 m n.m.
Střední sklon svahů
5,14°
- 21 -
3.2 Sběr dat Geodetická měření terénu probíhala ve dvou fázích. V první fázi – v září 2011 – prováděly dvě měřické skupiny zároveň měření dvěma různými metodami. První zvolenou metodou bylo tachymetrické měření pomocí totální stanice Topcon 105N v souřadném systému S-JTSK. Tato metoda se potýkala s několika problémy: špatná možnost připojení polygonového pořadu na body podrobného polohového bodového pole (nejbližší body PPBP se nachází na 2 km vzdáleném Blatném vrchu a na 2 km vzdálené Medvědí hoře), nevhodně zvolený polygonový pořad a zejména špatná viditelnost mezi odrůstajícími mladými smrky a suchými stojícími kmeny, kvůli které nebylo možno zachycovat výraznější lomové hrany a měřením pokrýt celé povodí. Proto byla tato metoda vyhodnocena jako nevhodná a pro potřeby diplomové práce nebyla takto získaná neúplná sada dat využita. Paralelně bylo prováděno druhé měření metodou diferenční GPS za použití GPS stanice Leica 1200. Pro korekci byla použita aktuální data z permanentních referenčních stanic sítě CZEPOS získávaná on-line přenosem mobilní internetovou sítí. Vzhledem k odlehlé poloze lokality u hranic s Německem, nebylo celé území dostatečně pokryto mobilním signálem, a proto nebylo možné celé povodí dokonale zmapovat. Navíc na severozápadním okraji povodí byl vlivem velice malého sklonu svahu nesprávně proveden prvotní odhad průběhu rozvodnice v terénu a v důsledku toho nebylo zaměřeno kompletně celé povodí. V září 2012 proběhla druhá fáze měření, za účelem doplnění chybějících dat. Bylo navázáno na diferenční GPS měření. Tentokrát však byla ke korekci použita mobilní referenční stanice – druhý přístroj Leica umístěný nad bodem o známých souřadnicích. Tím byl vyřešen problém s nedostatečným mobilním porytím a mohlo být dokončeno měření i na německé straně lokality. Celkově bylo během obou fází naměřeno pomocí georeferenční GPS zaměřeno 410 bodů, z čehož 38 tvoří nezávislý verifikační soubor.
- 22 -
3.3 Příprava dat Data naměřená pomocí GPS stanice byla extrahována pomocí softwaru Leica dodávaného s přístrojem. Pomocí stejného softwaru byla také data transponována z výchozího souřadnicového systému WGS 84 do systému S-JTSK. Na základě souřadnicového systému S-JTSK byl vytvořen pro větší přehlednost relativní souřadný systém s počátečním bodem v souřadnicích Y -830 000 a X -1150000, souřadnice Z zůstává nezměněna. Dále bylo potřeba data upravit do vhodného formátu, výstupní formát byl pro programování nevhodný. K tomu posloužily programy Microsoft Office Excel 2003 a 2010, ve kterých byla data manuálně setříděna tabulky o šesti sloupcích (číslo bodu, souřadnice Y, souřadnice X, nadmořská výška Bpv, 3D přesnost deklarovaná měřicím přístrojem, poznámka) a následně byla manuálně odstraněna chybná data, která svými hodnotami výrazně nezapadala do normálu. Poslední úpravou bylo doplnění poznámky u dat, u kterých nebyla vytvořena již při terénním měření. Pole poznámky totiž nemůže zůstat nevyplněné. V tomto tvaru byla tabulka exportována do formátu .txt (s tabulátorem, coby oddělovačem sloupců a tečkou – oddělovačem desetinných míst) a uložena jako Data.txt do pracovního adresáře. Slouží jako vstupní soubor do programu.
- 23 -
3.4 Struktura programu Program se dělí na dvě základní části: tvorba modelu a analýza modelu. Obě tyto části se skládají z většího množství různě složitých funkcí, základní kostra je však velice jednoduchá a přímočará. Vše je znázorněno v diagramu na obrázku 8:
Obrázek 8: Schéma programu.
- 24 -
V počátku procesu probíhá příprava výpočetního prostředí. Nejprve jsou nainstalovány dva potřebné externí balíčky. Balíček rgl slouží pro zobrazování 3D grafiky a balíček tripack pro tvorbu trojúhelníkové sítě z jednotlivých bodů. Následně probíhá načtení vlastních funkcí z pracovního adresáře. Dále začíná samotná práce s daty – nejprve načtením dat z pracovního adresáře ve formátu popsaném v kapitole 3.3. Následně jsou vyřazeny body, které nesplňují podmínky zadané jako argumenty (AccLimit, DistLimit) ve funkci TINmodel. První podmínka se vztahuje k přesnosti měření udávané přístrojem (Body zaměřené s menší přesností než je zadaný limit jsou ze seznamu odstraněny.), druhá k vzájemné vzdálenosti mezi jednotlivými body. (Pokud je vzdálenost mezi dvěma body menší než zadaný limit, jsou tyto považovány za zdvojený bod a jeden ze záznamů je odstraněn.) Následná tvorba trojúhelníkové sítě probíhá Delaunyho triangulací (případně Delaunyho triangulací s povinnými hranami – constraint Delaunay triangulation) pomocí externího balíčku tripack; je používán algoritmus, který popisuje Renka (1996). Další tři části (odstranění bezodtokých oblastí, vytvoření a analýza seznamu trojúhelníků a analýza hran) budou podrobněji rozebrány v kapitole 3.5. Výstupem z první části programu jsou dvě grafická okna, v prvním z nich je model znázorněn ve 3D ve druhém je znázorněn ve 2D s barevně rozlišenými směry odtoku jednotlivých trojúhelníků.
Analýza modelu začíná zjištěním rozvodnice povodí se zadaným závěrným profilem, následuje výpočet plochy povodí. Oběma algoritmům se budu rovněž věnovat v následující kapitole. U všech trojúhelníků je vypočtena jejich plocha pomocí obdobného algoritmu jako u výpočtu obsahu polygonu (3.5.6). U okrajových trojúhelníků povodí je navíc na základě lomových bodů rozvodnice a uzlových bodů modelu vypočteno, jaká jejich část se nachází uvnitř povodí. Na základě výsledků je dále váženým průměrem spočítán střední sklon svahů povodí a celková orientace terénu. Následuje verifikace modelu pomocí kritérií MAE a RMSE.
- 25 -
3.5 Vybrané algoritmy a výpočty Podrobným popisem každé funkce programu by rozsah této práce příliš nabyl na objemu, proto jsou v této kapitole vybrány a podrobněji popsány pouze některé důležité a zajímavé algoritmy, jejichž průběh není zřejmý na první pohled.
3.5.1 Odstranění bezodtokých oblastí Odstranění bezodtokých oblastí je bezpodmínečně nutné k vytvoření tzv. hydrologicky korektního modelu, tedy aby z každého bodu na povrchu modelu bylo možné vést křivku odpovídající trajektorii tekoucí povrchové vody až k okraji – konvexní obálce – modelu. Pro vyhledání bezodtokých oblastí je třeba nejprve nalézt lokální minima, tzv. dna bezodtokých oblastí. Pro určení dalších bodů bezodtokých oblastí je použito prohledávání do šířky. Označené bezodtoké oblasti jsou odstraněny změnou hodnoty souřadnice Z každého označeného bodu na hodnotu průměru hodnot Z jeho sousedních bodů. Pokud je bezodtoká oblast rozlehlejší, je nutné tento postup opakovat ve více iteracích. Algoritmus tedy vypadá následovně:
1. Vyhledej buňky, které mají menší hodnotu Z než všechny jejich sousední buňky a zároveň neleží na konvexní obálce trojúhelníkové sítě. Ulož je jako buňky náležící k bezodtokým oblastem (Drainless) a přidej je také do množiny měněných bodů (ChangedPoints). 2. Jejich sousedy zařaď do fronty. 3. Odeber z fronty první záznam a prozkoumej, zda daná buňka má kromě bodů obsažených v seznamu Drainless i jiné sousedy s nižší hodnotou souřadnice Z. Pokud nemá, zařaď jí do seznamu Drainless, přidej do množiny ChangedPoints a její sousedy, které nejsou obsaženy v seznamu Drainless zařaď na konec fronty. 4. Krok 3. opakuj, dokud není fronta prázdná. 5. Bodům označeným v seznamu Drainless změň souřadnici Z na hodnotu vypočtenou jako průměr souřadnic Z sousedů každého z nich. 6. Proveď krok 1. a pokud není nový seznam Drainless prázdný opakuj i body 2. až 6.
- 26 -
Zmíněný algoritmus je však vhodný pouze k odstranění menších bezodtokých oblastí uvnitř modelu. Při okrajích modelu (zejména v místě, kde odtok z povodí dosahuje okraje modelu) vytváří vady triangulace hráze, které neodpovídají realitě v terénu. Bezodtoké oblasti vzniklé těmito hrázemi je nutné řešit manuálně přidáním jednoho či několika bodů – tzv. technických bodů – vhodně umístěných do modelu přibližně v místech, kde očekáváme údolnici. Souřadnici Z u technických bodů volíme tak, aby voda z modelu mohla volně odtékat.
3.5.2 Vytvoření a analýza seznamu trojúhelníků Seznam trojúhelníků je vytvořen pomocí funkce triangles externího balíčku tripack. Každý trojúhelník je poté analyzován, je zjištěn jeho normálový vektor, na jehož základě je následně určen sklon trojúhelníku. Sklon trojúhelníku φ odpovídá odchylce roviny trojúhelníku od vodorovné roviny. Tato odchylka je rovna odchylce normálových vektorů zmíněných rovin. Sklon je tedy počítán pomocí rovnice odchylek dvou vektorů v prostoru:
cos ϕ =
a1b1 + a 2 b2 + a3b3 a1 + a 2 + a 3 2
2
2
b1 + b2 + b3 2
2
2
.
(5)
Pokud do rovnice dosadíme za a souřadnice normovaného normálového vektoru a za b normálový vektor vodorovné roviny b = (0,0,1), můžeme vhledem k tomu, že se jedná o dva jednotkové vektory rovnici zjednodušit na: cos ϕ = a3 .
(6)
Na základě normálového vektoru je vypočtena i směrová orientace trojúhelníku. Tentokrát je počítána odchylka průmětu normálového vektoru trojúhelníku do vodorovné roviny od jednotkového vektoru (0,1,0). Rovnice odchylky vektorů nelze zjednodušit jako v předchozím případě. Navíc v případě, že vektor promítnutý do vodorovné roviny směřuje do II. nebo III. kvadrantu, je potřeba vypočtenou hodnotu odečíst od plného úhlu.
- 27 -
Algoritmus analýzy je následující:
1. Vytvoř vektor U jako rozdíl 2. a 1. bodu trojúhelníku. 2. Vytvoř vektor V jako rozdíl 3. a 1. bodu trojúhelníku. 3. Vypočti souřadnice normálového vektoru jako vektorový součin vektorů U aV podle rovnice: u × v = (u 2 v3 − v 2 u 3 , u 3 v1 − v3u1 , u1v 2 − v1u 2 ) . 4. Pokud je hodnota Z normálového vektoru záporná, vynásob celý vektor -1. 5. Normálový vektor normuj. 6. Zjisti odchylku normálového vektoru trojúhelníku od normálového vektoru vodorovné roviny (0,0,1), převeď z radiánů na stupně. 7. Zjisti odchylku průmětu normálového vektoru trojúhelníku do vodorovné roviny od vektoru (0,1,0), převeď z radiánů na stupně. 8. Pokud souřadnice X normálového vektoru je menší než nula, odchylku vypočtenou v kroku 7 odečti od 360. 9. Opakuj body 1 až 8 pro každý trojúhelník.
- 28 -
3.5.3 Analýza hran V modelu rozlišuji čtyři typy hran trojúhelníků podle sklonu trojúhelníků, jimž náleží a které oddělují. Prvním typem jsou hrany, které tvoří konvexní obálku. Každá z těchto hran náleží pouze jednomu trojúhelníku a je proto snadné je odlišit od ostatních typů. Všechny ostatní hrany náleží dvěma trojúhelníkům zároveň a v závislosti na sklonu trojúhelníků je můžeme rozdělit na hřbetnice, údolnice a hrany ve svahu. Toto rozdělení a označení jednotlivých hran je velice důležité pro následné analýzy tras odtoku povrchové vody a pro hledání rozvodnice. Algoritmus funguje na porovnání Z souřadnic tří význačných bodů ležících na kolmici ke zkoumané hraně a vypadá takto:
1. Prozkoumej sousední trojúhelníky, pokud existuje jen jeden sousední trojúhelník, přiřaď hraně typ 0 (konvexní obálka), kroky 2 až 4 přeskoč a pokračuj krokem 5. 2. V průmětu do vodorovné roviny najdi střed zkoumané hrany. Tímto bodem veď přímku kolmou ke zkoumané hraně a najdi její průsečíky se zbylými hranami obou sousedních trojúhelníků. 3. Tři význačné body nalezené v kroku 2 (střed zkoumané hrany a dva průsečíky) promítni zpět na povrch modelu, zjisti Z souřadnice každého z nových bodů. 4. Zjištěné souřadnice porovnej: • Pokud Z souřadnice středu zkoumané hrany je vyšší než souřadnice obou průsečíků, přiřaď hraně typ 1 (hřbetnice) • Pokud Z souřadnice středu zkoumané hrany je nižší než souřadnice obou průsečíků, přiřaď hraně typ 2 (údolnice) • Pokud se Z souřadnice středu zkoumané hrany nachází mezi Z souřadnicemi průsečíků, přiřaď hraně typ 3 (hrana ve svahu) 5. Kroky 1 až 4 proveď pro každou hranu.
- 29 -
Na obrázku Obrázek 9 jsou znázorněny popisované situace ve volném rovnoběžném zobrazení se zvýrazněnými význačnými body a jejich průměty do roviny π:
Obrázek 9: Typy hran: a) hřbetnice, b) údolnice a c) hrana ve svahu.
- 30 -
3.5.4 Nalezení trasy odtoku Algoritmus nalezení trasy odtoku je jeden z nejdůležitějších v celém programu. Je mnohokrát využíván i v dalších částech programu. Ze zadaných parametrů - souřadnic Y a X - nalezne nejstrmější cestu na okraj modelu, tato cesta simuluje trajektorii odtoku povrchových vod. Celý výpočet probíhá v průmětu do vodorovné roviny a pro jednotlivé lomové body trajektorie je dopočtena souřadnice Z odpovídajícího bodu na povrchu terénu. Výpočet celé trajektorie začíná v zadaných souřadnicích a probíhá v cyklech, přičemž v každém cyklu se posune do nového lomového bodu trajektorie (tzv. aktuální pozice) přes jeden trojúhelník ve směru spádnice na některou z protějších hran, nebo po jedné hraně do nižšího uzlového bodu modelu. 1. Zjisti číslo trojúhelníku, ve kterém se nachází zadané souřadnice. 2. V průmětu do vodorovné roviny nalezni průsečíky přímky dané zadaným bodem a normálovým vektorem trojúhelníku s hranami tohoto trojúhelníku. 3. Vyber ten průsečík, který odpovídá bodu s nižší souřadnicí Z a ulož ho jako tzv. aktuální pozici. 4. Zaznamenej současnou aktuální pozici jako poslední bod trajektorie. 5. Pokud aktuální pozice odpovídá některému z uzlových bodů modelu, pokračuj krokem 6a, pokud neodpovídá žádnému uzlovému bodu, pokračuj krokem 6b. 6. a) Nalezni všechny směry odtoku z aktuální pozice – údolnice a spádnice trojúhelníků – a vyber tu s největším kladným sklonem. Pokud se jedná o hranu – údolnici, pokračuj krokem 7a, pokud se jedná o spádnici trojúhelníku pokračuj krokem 7b. b) Pokud je typ hrany (viz. kapitola 3.5.3), na které se nachází aktuální pozice, údolnice, pak pokračuj krokem 7a, jinak pokračuj krokem 7b. 7. a) Přepiš aktuální pozici na bod určený souřadnicemi nižšího bodu hrany. Pokračuj krokem 8. b) Urči číslo trojúhelníku, přes který poteče voda, a přepiš aktuální pozici na nižší z průsečíků hran tohoto trojúhelníku a přímky určené vektorem jeho spádnice a dosavadními souřadnicemi aktuální pozice. 8. Pokud se aktuální pozice nenachází na konvexní obálce modelu, opakuj kroky 4 až 8.
- 31 -
3.5.5 Nalezení rozvodnice Algoritmus hledání rozvodnice funguje na velmi podobném principu jako algoritmus hledání trasy odtoku. Vstupním parametrem je však tentokrát řetězec, který se shoduje s poznámkou ve vstupním souboru dat u bodu odpovídajícího uzávěrovému profilu hledaného povodí. Algoritmus probíhá ve dvou větvích, které začínají v místě uzávěrového profilu a opět se stýkají v nejvyšším bodě povodí.
1. Zjisti číslo bodu, u kterého se shoduje poznámka s řetězcem vstupního parametru funkce. 2. Najdi lokální maxima modelu, body jejichž všichni sousedé mají nižší hodnotu souřadnice Z. 3. Najdi trojúhelníky obsahující uzlový bod zjištěný v kroku 1. Z těžiště každého z nich zjisti trasu odtoku. Rozděl trojúhelníky na takové, z jejichž těžiště vede trasa odtoku přes uzávěrový profil, a na ty, ze kterých nevede. 4. Najdi dvě tzv. kritické hrany – hrany na rozmezí mezi trojúhelníky, jejichž těžiště v jednom případě náleží do povodí, ve druhém ne. 5. Pro každou kritickou hranu zjisti, zda je její typ hřbetnice. Pokud ano, jako výchozí bod dané větve rozvodnice urči okrajový bod hrany s vyšší hodnotou souřadnice Z. Pokud kritická hrana není hřbetnicí, najdi průsečíky hran obou sousedních trojúhelníků s přímkou danou souřadnicemi uzávěrového profilu a spádnicí každého z trojúhelníků. Jako výchozí bod větve rozvodnice v tomto případě vyber jediný průsečík ze čtyř, který se neshoduje s uzávěrovým profilem. 6. Výchozí bod větve ulož jako aktuální pozici. 7. Aktuální pozici ulož jako poslední lomový bod rozvodnice. 8. Pokud je aktuální pozice shodná se souřadnicemi některého z uzlových bodů modelu pokračuj krokem 9a, jinak pokračuj krokem 9b 9. a) Obdobným způsobem jako v bodech 3 a 4 najdi kritické hrany. Vyber tu, která směřuje vzhůru, a ověř, zda je hřbetnicí. Pokud je hřbetnicí, pokračuj krokem 10a, pokud není, pokračuj krokem 10b. b) Pokud je hrana, na které se nachází aktuální pozice hřbetnicí, pokračuj krokem
- 32 -
10a, pokud není, pokračuj krokem 10b. 10. a) Jako novou aktuální pozici ulož souřadnice okrajového toho bodu s vyšší hodnotou souřadnice Z hrany, na které se nachází dosavadní pozice (příp. kritické hrany). Pokračuj krokem 11. b) Jako novou aktuální pozici ulož průsečík přímky dané aktuální pozicí a vektorem spádnice sousedního trojúhelníku s jeho protější hranou. 11. Pokud se aktuální pozice neshoduje s některým z lokálních maxim, opakuj kroky 8 až 10. 12. Kroky 6 až 11 opakuj pro obě větve rozvodnice. 13. Polygon rozvodnice vytvoř jako vektor bodů, který začíná uzávěrovým profilem, pokračuje body první větve rozvodnice a body druhé větve rozvodnice zkrácené o poslední bod – lokální maximum, v převráceném pořadí.
3.5.6 Výpočet plochy povodí Plocha povodí je vypočtena z polygonu rozvodnice pomocí rovnice (Beyer 1987):
x1 y1 A = abs
x2 y2
+
x2
x3
y2
y3
+ ... +
2
xn yn
x1 y1 ,
(7)
kde n je počet lomových bodů rozvodnice. Pro potřebu programu byl vzorec upraven následovně:
x n x1 xi xi +1 y n y1 n −1 y i yi +1 A = abs +∑ 2 2 i =1
.
(8)
Vlastní algoritmus lze z této úpravy vzorce snadno vyvodit.
- 33 -
3.5.7 Verifikace modelu Model je verifikován pomocí nezávislé verifikační sady měřených bodů. Pro každý verifikační bod je zjištěna modelovaná nadmořská výška ve stejných souřadnicích Y a X. Na takto získaná rezidua je pohlíženo skrze kritéria MAE a RMSE. Vzorce kritérií jsou dány rovnicemi:
n
MAE =
∑ i =1
Z mod [i ] − Z ver [i ] (9)
n
a
∑ (Z n
RMSE =
i =1
mod [i ]
− Z ver [i ]
n
)
2
,
(10)
ve kterých Zmod značí souřadnici Z určenou modelem, Zver značí měřenou souřadnici Z ze souboru verifikačních dat a n je počet verifikačních bodů.
- 34 -
3.6 Práce s programem V této kapitole je krok po kroku uvedeno, jak pracovat s programem přiloženým na CD. Příkazy zadávané do konzole programu R jsou v textu uvedeny na jednotlivých řádcích a jsou předznamenány symbolem „>“. Tento symbol se do příkazového řádku nezadává. Nezadává se ani případné interpunkční znaménko uvedené na konci řádku s příkazem, příkaz musí vždy končit symbolem „)“
3.6.1 Příprava prostředí Pro práci s programem je třeba mít nainstalovaný software R (volně dostupný ze stránek projektu www.r-project.org). Všechny skripty byly napsány pro verzi R 2.12.0, pro jiné verze nelze funkčnost skriptů zcela garantovat. Před započetím práce je potřeba připravit pracovní adresář, do kterého je nutné nahrát všechny skripty (soubory formátu *.R) z adresáře Skripty na přiloženém CD a vstupní soubor dat ve formátu *.txt v potřebné úpravě (viz kapitola 3.3). Pokud máme soubory připravené a software nainstalovaný, spustíme jej poklepáním na jeho ikonu. Objeví se okno s konzolou R Console, ve kterém je potřeba nastavit cestu k pracovnímu adresáři příkazem:
>setwd(“adresa”),
kam je navíc zapotřebí místo slova adresa zkopírovat cestu k připravenému pracovnímu adresáři. Příkaz vepíšeme či zkopírujeme do posledního, aktivního řádku konzole vyplníme adresu a potvrdíme klávesou enter. Nastavení pracovního adresáře lze rovněž provést v rozbalovacím menu File – Change dir, kde správný adresář vybereme ze zobrazeného adresářového stromu. Před započetím práce je také vhodné odstranit všechny proměnné, které mohou být v paměti programu uloženy z předchozích sezení. K tomu slouží příkaz:
>rm(list=ls()).
- 35 -
Pokud již v počítači nainstalované nejsou, je dále nutné instalovat dva externí balíčky používané v předkládaném programu. V rozbalovacím menu klikneme na Packages – Instal packages a v okně Packages vybrat možnosti rgl a tripack. Vícenásobný výběr se provádí se stisknutou klávesou ctrl. Je možné, že před samotným výběrem balíčku, bude nutné vybrat také umístění serveru, odkud mají být balíčky staženy (tzv. CRAN mirror). V takovém případě se namísto okna Packages nejdříve zobrazí okno Select CRAN mirror, ve kterém provedeme výběr serveru. Nejvhodnější je zvolit umístění serveru v některé blízké lokalitě. Poslední přípravnou akcí je načtení všech funkcí programu. To provedeme zadáním příkazu:
>source(“sourceAll.R”).
3.6.2 Tvorba modelu terénu Model terénu vytvoříme ze zadaného souboru dat vyvoláním funkce TINmodel:
>Model=TINmodel(),
čímž zároveň do proměnné Model (typ list) uložíme všechny potřebné informace o vytvořeném modelu terénu. Název proměnné můžeme samozřejmě vybrat libovolný jiný. Do kulatých závorek je možné zadávat hodnoty následujících argumentů oddělené čárkami. DataSource – udává název vstupního souboru dat v pracovním adresáři. Zadává se ve formátu DataSource=“Název.txt“. AccLimit – udává limit požadované přesnosti měření deklarované stanicí GPS v metrech. Zadáním hodnoty 0 nebo záporné hodnoty bude třídění vynecháno a budou použita všechna data bez ohledu na přesnost měření. DistLimit – udává limit minimální vzájemné vzdálenosti dvou uzlových bodů v metrech. Body, které jsou k sobě blíže, budou považovány za jeden zdvojený bod a záznam s menší přesností měření bude odstraněn. Zadáním hodnoty 0 nebo záporné hodnoty bude třídění vynecháno. ConstraintNote – je vektor proměnných typu string, které označují uzlové body, jež jsou určeny jako uzlové body tzv. povinných hran Delaunayho triangulace. Tyto
- 36 -
body musí být ve vstupním souboru dat uvedeny ve správném pořadí a se shodným označením
v posledním
sloupci
Poznámka.
Zadává
se
ve
formátu
ConstraintNote=c(“Poznámka 1”,”Poznámka 2”,”Poznámka 3”). Pokud je vektor zanechán prázdný, v modelu nebudou použity povinné hrany. Pokud jsou některé poznámky
v příkazu
vynechány,
nabývají
těchto
výchozích
hodnot:
DataSource="Data.txt", AccLimit=0.15, DistLimit=0.15, ConstraintNote=c(). Výpočet může při vyšším počtu záznamů vstupního souboru trvat i několik minut. Po úspěšně provedeném výpočtu se objeví nové okno grafického rozhraní rgl s trojrozměrným znázorněním modelu terénu a grafické okno s dvourozměrným znázorněním modelu terénu s barevně odlišenými trojúhelníky různé orientace, doplněné o okno s legendou. Posledním krokem tvorby modelu je extrahování dat z proměnné Model a jejich uložení jako globální proměnné:
>globalize(Model).
Tím je model hotov a připraven k analýze.
3.6.3 Analýza modelu Funkci pro analýzu modelu vyvoláme příkazem:
>Basin=basinAnalysis(Model,BEnote),
kde do kulatých závorek uvedeme až dva parametry. Parametr Model je povinný, zadáváme jím proměnnou, do které byl uložen výstup z procedury TINmodel (v našem případě proměnná Model), a na jejímž základě je prováděna analýza. Druhým parametrem je BEnote, který udává řetězec (ohraničený uvozovkami), jímž je označen uzávěrový profil povodí v poznámce vstupního souboru dat. Pokud druhý parametr vynecháme, bude použita výchozí hodnota BEnote=“PRELIV“. Procedura je poněkud náročná, a proto může výpočet trvat delší dobu. Po výpočtu je do aktivního grafického okna vykreslena rozvodnice. Informace o povodí – výstup z analýzy – lze zobrazit v konzoli vypsáním proměnné BasinInfo.
- 37 -
>BasinInfo
Proměnná BasinInfo obsahuje informace o souřadnicích uzávěrového profilu, ploše, průměrném sklonu a orientaci povodí.
3.6.4 Interaktivní funkce V programu je obsaženo i několik interaktivních funkcí. Ty se vždy aktivují zadáním příslušného příkazu v příkazovém řádku konzole. Tedy:
>zoomIN() >thisPoint() >nearestNode() >placeDrop() nebo >distance().
Následně se spustí po kliknutí do příslušného místa v aktivním grafickém okně.
Funkce zoomIN otevře nové grafické okno a vytiskne do něj výřez modelu určený dvěma diagonálně umístěnými rohovými body načtenými dvojím kliknutím do původního grafického okna. Nově otevřené grafické okno se okamžitě stává aktivním, a proto všechny další interaktivní funkce se vztahují právě k němu. Zavřením okna s výřezem se stane opět aktivním předcházející grafické okno. Do kulatých závorek je možné uvést parametr type, který specifikuje typ zobrazení modelu. Type=1 barevným odlišením trojúhelníků vykresluje orientaci terénu, Type=3 vykresluje sklon jednotlivých trojúhelníků, Type=2 je kombinací předchozích možností.
Funkce thisPoint zobrazí informace o bodu označeném kliknutím do zobrazení modelu. Zobrazovanými informacemi jsou Y, X a Z souřadnice bodu, číslo trojúhelníku, kterému zadaný bod náleží, a informace o tomto trojúhelníku (sklon, orientaci terénu,
- 38 -
plochu trojúhelníku – resp. jeho průmětu do vodorovné roviny a plochu, kterou daný trojúhelník přispívá do povodí).
Funkce nearestNode zobrazí informace o nejbližším uzlovém bodu modelu (číslo bodu v modelu, Y, X a Z souřadnice a poznámka k bodu uvedená ve zdrojovém souboru dat).
Funkce distance vypíše do konzole vzdálenost (v metrech) dvou bodů označených v modelu.
Poslední interaktivní funkcí je placeDrop, která zobrazí trasu povrchového odtoku z příslušného bodu modelu zadaného kliknutím do aktivního grafického okna.
3.6.5 Uložení a načtení modelu Program umožňuje hotový model uložit a následně znovu načíst a zobrazit, což je výhodné vzhledem k náročnosti výpočtu. K ukládání projektu slouží funkce saveP:
>saveP(File,ToSave),
kde se do kulatých závorek uvádí dva parametry. Prvním je File, název souboru, do kterého bude projekt uložen. Řetězec názvu se zadává v uvozovkách. Pokud je tento parametr vynechán, je použita defaultní hodnota File=“project1“. Druhým parametrem je seznam proměnných, které mají být uloženy,
ToSave. Zadává se ve formátu
ToSave=list(Model,Basin). K načtení projektu slouží funkce
>loadP(File),
kde se za parametr File zadává název načítaného souboru v uvozovkách. Po načtení souboru jsou otevřena grafická okna se znázorněním modelu ve 2D a v 3D a program je připraven k dalšímu užívání.
- 39 -
3.6.6 Funkce pro ovládání grafiky Funkce pro grafické zobrazení modelu jsou integrovány už i v obalových funkcích TINmodel a basinAnalysis, nicméně je možné je používat i samostatně. Je však možné je použít pouze v případě, že základní model byl již vytvořen. Dvourozměrný model je možné zobrazit funkcí:
>plotModel2D(Points,Triangles,TSlope,TDirection,Type=1,MaxSlope=30).
Parametry Points, Triangles, TSlope a TDirection necháváme vyplněny jako v uvedeném příkladu, Type udává typ zobrazení (v závislosti na zobrazovaných veličinách) a může nabývat hodnot od 1 do 3 (1 = orientace terénu, 2 = orientace a sklon terénu, 3 = sklon terénu). Pokud zůstane parametr nezadán, je použita výchozí hodnota 1. Parametr MaxSlope má význam pouze při hodnotě předchozího parametru Type=3, ovlivňuje rozsah sklonů odlišovaných v modelu rozdílnými odstíny. Maximální sklon se zadává ve stupních, a pokud zůstane parametr nevyplněn, je použita výchozí hodnota 30°. Trojúhelníky se sklonem vyšším, než je udaný tímto parametrem, jsou vyplněny černou barvou a jejich sklony tedy nelze odlišit. Prostorové zobrazení modelu je spouštěno příkazem:
>plotModel3D (Points,Triangles,ChangedPoints).
Otevírá se v samostatném okně rgl a je možné jej ovládat myší (přiblížení a rotace modelu). Funkce zobrazení rozvodnice do modelu může proběhnout pouze v případě, že je průběh rozvodnice již spočítán. Je vyvolávána příkazem:
>plotBasinDivide().
Poslední funkcí pro ovládání grafiky je:
>plotRidgeEtThalweg(Points,Arcs,Triangles,TSlope,TDirection).
- 40 -
Zobrazuje do dvourozměrného modelu, která hrana je údolím a která hřebenem. Hřebeny jsou zvýrazněny silnější černou linkou, údolí šedě. Funkce není obsažena v balíku načítaném příkazem sourceAll (viz Příprava prostředí 3.6.1), a proto je nutné ji před případným použitím načíst příkazem:
>source(“plotRidgeEtThalweg.R”).
3.6.7 Příklad použití Na následujících řádcích je na měřených datech z povodí Modrava1 uveden příklad, jak může vypadat sekvence příkazů zadávaná do příkazového řádku konzole.
>Model=TINmodel (DataSource=”Data.txt”, AccLimit=0.15, DistLimit=0.15) >globalize(Model) >Basin=basinAnalysis(Model,BEnote=”PRELIV”) >BasinInfo >saveP(File=”Priklad”, ToSave=list(Model,Basin))
V příkladu je nejdříve vytvořen celkový model, který je uložen do proměnné Model. Všechny jeho důležité části jsou dále uloženy do globálních proměnných, na jejichž základě je provedena analýza povodí, jehož uzávěrový profil je označen poznámkou „PRELIV“. Dále jsou v konzoli vytištěny informace o povodí a celý projekt je uložen do souboru s názvem Priklad. Při dalším spuštění programu R je možné tento projekt opět načíst příkazem:
>loadP(“Priklad”).
- 41 -
4 Výsledky a výstupy digitálního modelu terénu Kapitoly Výstupy analýzy terénu i Výstupy analýzy povodí předkládají výsledky pouze z původního, nezredukovaného souboru dat Model M1, výstupy z ostatních modelů jsou vyobrazeny v přílohách 2 - 5. V kapitole Výsledky verifikace a porovnání modelů jsou modely porovnány z hlediska kritérií MAE (mean absolut error - střední absolutní chyba) a RMSE (root mean square error – odmocněná střední kvadratická chyba) (viz rovnice 7 a 8). Výstupem z digitálního modelu terénu jsou převážně obsáhlé datové tabulky, které nejsou v textové podobě vůbec vypovídající, proto jsou pro prezentaci výsledků použita zejména grafická znázornění pomocí obrázků, grafů a diagramů.
- 42 -
4.1 Výstupy analýzy terénu 4.1.1 Bodové pole V terénu bylo celkem zaměřeno 410 bodů, z toho 38 tvoří samostatný soubor určený pro verifikaci modelu. Část základního souboru (7 bodů) byla odstraněna již při manuální přípravě dat. Ze souboru dat určených pro tvorbu modelu bylo 72 položek odstraněno pro nedostatečnou přesnost měření deklarovanou přístrojem (při limitu přesnosti 15 cm). Další dva body byly posléze odstraněny pro jejich duplicitu s jinými. (Limitní vzdálenost dvou bodů byla stanovena rovněž na 15 cm.) Tři body byly naopak do souboru manuálně přidány při odstraňování bezodtokých oblastí. Model je tedy tvořen 294 uzlovými body, které jsou zhuštěny v místech s větší členitostí terénu a v místech, kde byl předpokládán průběh rozvodnice a údolnice zkoumaného povodí. Jejich rozmístění je zobrazeno na obrázku.
Obrázek 10: Bodové pole modelu M1.
- 43 -
4.1.2 Trojúhelníková síť Trojúhelníková síť byla vytvořena Delaunyho triangulací bez povinných hran. Celý model je tvořen 571 trojúhelníky s celkem 864 hranami. Konvexní obálku modelu tvoří 15 bodů. Celková plocha modelu čítá 38,57 ha.
Obrázek 11: Triangulační síť modelu M1.
- 44 -
4.1.3 Odstranění bezodtokých oblastí Vlivem vad triangulace při okraji modelu vznikla v severní části modelu (v místě kde se tok opouští model) hráz, která neodpovídá reálnému terénu a díky níž vznikla v modelu rozsáhlá bezodtoká oblast. Pokud by taková oblast byla odstraněna prostým vyplněním pomocí algoritmu uvedeného v kapitole 3.5.1, došlo by ke značnému znehodnocení modelu. Proto byly do zdrojových dat manuálně přidány tzv. technické body, které tuto „hráz“ odbouraly, a model byl vytvořen znovu z aktualizovaných vstupních dat. V těchto místech proto model sice neodpovídá realitě, ale vzhledem k tomu, že se tato část nachází až pod uzávěrovým profilem zájmového povodí, nebude výpočet rozvodnice ovlivněn. Hráz je patrná jako červené trojúhelníky v severní části modelu na obrázku v příloze 6, obrázek znázorňuje orientaci svahů modelu před doplněním technických bodů. Souřadnice vkládaných technických bodů jsou uvedeny v tabulce 2.
Tabulka 2: Souřadnice technických bodů. Pořadové číslo bodu 292 293 294
Označení bodu 9997 9998 9999
Souřadnice Y (East) -1430 -1415 -1410
Souřadnice Nadmořská X (North) výška (Bpv)[m] -5080 1212 -5070 1211 -5050 1210
Na nový model byl aplikován výše zmíněný algoritmus pro automatické odstranění bezodtokých oblastí, žádná další bezodtoká oblast však nebyla nalezena, a proto u žádného bodu nemusela být změněna souřadnice Z.
- 45 -
4.1.4 Orientace trojúhelníků Orientace trojúhelníků a směry odtoku pro každý trojúhelník byly vypočteny na základě metodiky uvedené v kapitole 3.5.2, výsledek je graficky prezentován na následujícím obrázku. Pro rozlišení orientace trojúhelníku byla použita paleta 360 odstínů, každý odstín pro jednostupňový interval. Jednotlivé odstíny a směry, které symbolizují, jsou zobrazeny v kruhovém diagramu dole na obrázku 12.
Obrázek 12: Model M1 - orientace terénu a směry odtoku.
- 46 -
V grafu 1 je znázorněno rozložení orientace terénu v celém modelu v desetistupňových intervalech.
Graf 1: Rozdělení orientací terénu v modelu M1.
Z obou znázornění orientací je patrné, že model zobrazuje hlavně svah se severovýchodní až východní orientací. Data rozložení orientací terénu modelu M1 v tabulkové podobě jsou součástí přílohy 2.
- 47 -
4.1.5 Sklony trojúhelníků Sklon terénu v jednotlivých trojúhelnících byl vypočten metodou uvedenou v kapitole 3.5.2. Grafické znázornění sklonů mezi 0 a 30 stupni je zobrazeno na obrázku 13.
Sklon [°]
Obrázek 13: Model M1 - sklony svahů.
Průměrný sklon svahů v celém modelu je 6,17°. Terén je v lokalitě velmi málo členitý, výraznějším prvkem je pouze pramenná jímka v severní části. Trojúhelníky v blízkosti okraje modelu sklonem nemusí odpovídat skutečnému terénu vzhledem k jejich
- 48 -
neideálnímu, protáhlému tvaru. Rozložení sklonů je znázorněno v grafech 2 a 3. První z grafů vykresluje rozložení sklonů v normálním měřítku, druhý s logaritmickým měřítkem na ose y. Zdrojové tabulky grafů jsou součástí přílohy 2.
Graf 2: Rozdělení sklonů terénu v modelu M1v normálním měřítku.
Graf 3: Rozdělení sklonů terénu v modelu M1 v logaritmickém měřítku na ose y.
- 49 -
4.1.6 Analýza hran Analýza hran je prováděna zejména pro potřeby dalšího průběhu programu, samotná není příliš zajímavým ani vypovídajícím výstupem. Výsledky jsou znázorněny na obrázku v příloze 2.
4.2 Výstupy analýzy povodí Předmětem analýzy modelu byl terén experimentálního povodí Modrava 1, s uzávěrovým profilem situovaným v místě hrotu thomsnova měrného přelivu. Souřadnice tohoto bodu v systému S-JTSK jsou Y = -831450.068 a X = -1155107.998. Nadmořská výška v systému Balt po vyrovnání je 1212.642 m. Nejprve byl zjištěn tvar povodí, z něj celková plocha povodí a dalšími výpočty pak tvarové charakteristiky povodí. V této kapitole uvádím výsledky modelu M1, který je brán jako stěžejní.
- 50 -
4.2.1 Rozvodnice Hledání rozvodnice probíhalo podle algoritmu uvedeného v kapitole 3.5.5 – od uzávěrového profilu cestou největších sklonů, tj. nejstrmějších spádnic trojúhelníků, a hřebenů, až k vrcholu povodí. Práce byla usnadněna specifickým tvarem povodí, které je tvořeno pouze jedním svahem, rozvodnice proto probíhá pouze přes jeden vrchol – lokální maximum – a neprobíhá žádným sedlem. Průběh rozvodnice je patrný z obrázku 14.
Obrázek 14: Model M1 – rozvodnice.
- 51 -
4.2.2 Tvarové charakteristiky povodí Na základě vzorce z kapitoly 3.5.6 byla určena plocha povodí, následně byl zjištěn průměrný sklon svahů a celková orientace povodí. Výsledné hodnoty jsou uvedeny v tabulce 3.
Tabulka 3: Hodnoty terénních charakteristik povodí v modelu M1.
Plocha povodí [m2]
102919
Délka rozvodnice [m]
1327.1
Délka povodí [m]
580.9
Průměrný sklon svahů [°]
6.68
Celková orientace povodí [°]
64.99
Min. nadmořská výška [m n.m.]
1212,6
Max. nadmořská výška [m n.m.]
1279,2
Orientace svahů povodí Povodí sestává prakticky pouze z jediného svahu, tomuto faktu odpovídá i rozložení orientací svahů, kde naprostá většina plochy leží mezi 50 a 90 stupni. Výsledky jsou zobrazené v grafu 4 a v tabulce 4.
- 52 -
Graf 4: Rozdělení orientací terénu v povodí pro model M1.
Tabulka 4: Hodnoty rozdělení orientací terénu v povodí pro model M1. Orientace[°] 0-10 10-20 20-30 30-40 40-50 50-60 60-70 70-80 80-90 90-100 100-110 110-120 120-130 130-140 140-150 150-160 160-170 170-180
Plocha[m^2] 325 411 3430 5774 6547 8231 35177 26077 11329 4223 36 131 63 35 0 50 29 0
Orientace[°] 180-190 190-200 200-210 210-220 220-230 230-240 240-250 250-260 260-270 270-280 280-290 290-300 300-310 310-320 320-330 330-340 340-350 350-360
Plocha[m^2] 0 0 0 0 0 0 0 0 0 18 0 0 0 0 0 39 125 334
- 53 -
Sklon svahů povodí Terén povodí je jednotvárný, mírně svažitý. Jediný výraznější terénní prvek je pramenná jímka v severní části povodí. Většina plochy povodí je odkloněna od vodorovné roviny pod úhlem čtyř až deseti stupňů. Rozdělení sklonů je znázorněno v grafech 5 a 6 a v tabulce 5.
Graf 5: Rozdělení skolnů svahů v povodí pro model M1 v normálním měřítku.
- 54 -
Graf 6: Rozdělení sklonů svahů v povodí pro model M1 v logaritmickém měřítku na ose y.
Tabulka 5: Hodnoty rozdělení sklonů v povodí pro model M1. Sklon[°] 0-2 2-4 4-6 6-8 8-10 10-12 12-14 14-16 16-18 18-20 20-22 22-24 24-26 26-28 28-30 30-32 32-34 34-36 36-38 >38
Plocha[m^2] 2500 6995 30007 40569 19844 1543 120 247 118 34 224 14 112 0 1 26 0 0 30 0
- 55 -
4.3 Výsledky verifikace a porovnání modelů Každý z modelů byl verifikován sadou 38 nezávislých bodů podle kritérií MAE a RMSE, jejichž rovnice jsou uvedeny v kapitole 3.5.7. Cílem bylo určit, do jaké míry lze generalizovat povrch terénu v modelu povodí Modrava 1, aniž by došlo k výraznému poklesu pravdivosti modelu, a zda by naopak nebylo vhodné Model M1 doplnit o další měřené body. Hodnoty verifikačních kritérií jsou uvedeny v tabulce 6 a znázorněny v grafech 7 a) a 7 b).
Tabulka 6: Hodnoty verifikačních kriterií MAE a RMSE pro jednotlivé modely. M1 Počet uzlových bodů modelu MAE [m] RMSE [m]
Model M1_75 M1_50 M1_25 294 222 147 71 0.280 0.282 0.334 0.794 0.409 0.362 0.429 1.120
- 56 -
Graf 7: Hodnoty verifikačních kriterií: a) MAE, b) RMSE.
V obrázku 15 je zachyceno prostorové rozložení reziduí jednotlivých modelů.
- 57 -
Obrázek 15: Prostorové rozložení reziduí v jednotlivých modelech.
- 58 -
4.3.1 Porovnání modelů vytvořených na bázi trojúhelníkové sítě Čtyři vytvořené modely podávají čtyři různé sady výsledků. Pro jejich porovnání slouží tabulka nejdůležitějších charakteristik modelových povodí (Tabulka 1Tabulka 7).
Tabulka 7: Srovnání terénních charakteristik povodí udávaných jednotlivými modely.
Plocha povodí Délka rozvodnice Průměrný sklon svahů Orientace povodí
Model Jednotka M1 M1_75 M1_50 M1_25 2 [m ] 102919 94104 106876 81006 [m] 1327.07 1290.57 1347.95 1316.60 [°] 6.68 6.82 6.48 6.48 [°] 64.99 66.39 61.91 61.86
V obrázku 16 jsou vyznačeny, rozvodnice vytvořené jednotlivými modely. Nejzásadněji se od ostatních odlišuje rozvodnice vytvořená modelem M1_25, zbylé mají velice podobný průběh. Rozvodnice modelu M1_75 se poněkud výrazněji odlišuje pouze v horní části a na jihovýchodní hranici povodí, rozvodnice modelu M1_25 tvarově kopíruje základní model M1, pouze v některých místech se mírně odklání směrem ven z povodí.
- 59 -
Obrázek 16: Porovnání průběhu a tvaru rozvodnic nalezených různými modely.
4.3.2 Srovnávací model na bázi čtvercové sítě Model sloužící pro rámcové porovnání metod založených na bázi nepravidelné trojúhelníkové sítě a pravidelné čtvercové sítě, byl vytvořen v programu ArcGIS. Byl použit nejjednodušší algoritmus pro distribuci odtoku z jednotlivých bodů, algoritmus zvaný D8. Různé používané algoritmy shrnuje Bašta (2008). Plocha povodí vypočtená tímto modelem čítá 107 536 m2; další, grafické výstupy modelu na bázi čtvercové sítě jsou připojeny v příloze 7.
- 60 -
5 Diskuze Cílem této kapitoly je zhodnocení použitých metod uvedených v metodice v kapitole 3 a shrnutí a zhodnocení výsledků a výstupů diplomové práce uvedených v kapitole 4.
5.1 Měřické metody Pro sběr tachymetrických dat byly použity dvě techniky. První zvolenou technikou bylo použití totální stanice Topcon 105N, druhou použití referencované GPS stanice Leica 1200. Užití totální stanice se ukázalo být nevhodným, vzhledem ke značné míře pokrytí povrchu lokality odrůstajícím smrkovým podrostem. V podrostu bylo možné nacházet průhledy přes značnou část plochy povodí, nicméně zaměřovat konkrétní body, lomové hrany atp. potřebné pro model terénu na bázi nepravidelné trojúhelníkové sítě bylo téměř nemožné. Práce postupovala velmi pomalu, z vytvořeného polygonového pořadu navíc nebylo možné dobře zaměřit celou plochu povodí. Další nevýhodou byla velká vzdálenost bodů podrobného polohového bodového pole potřebných pro zaměření hlavního polygonového pořadu. Jako vhodnější se ukázalo použití georeferencované GPS stanice. Vzhledem k velmi odlehlé lokalitě a nedostatečnému pokrytí mobilním signálem však nebylo možné v celé ploše povodí použít informace ze stabilních referenčních stanic sítě CZEPOS. Proto byla v poslední fázi měření užita mobilní referenční stanice ustanovená nad bodem o známých souřadnicích. To se nakonec ukázalo jako nejlépe použitelná měřická metoda pro povodí Modrava 1. Výraznou výhodou této lokality pro tento způsob měření byla absence stromového patra a tím velice dobrá dohlednost družic GPS.
- 61 -
5.2 Digitální model terénu Cílem diplomové práce bylo vytvořit digitální model terénu povodí Modrava 1. Základním požadavkem na model je samozřejmě, aby co nejlépe reprezentoval skutečný terén, a zároveň si zachoval dostatečnou jednoduchost.
Generalizace terénu Ze základní sady vstupních dat byly vytvořeny čtyři modely, které využívaly různou poměrnou část vstupního souboru dat (100 %, 75 %, 50 % a 25 % bodů). Z hodnot vypočtených kritérií uvedených v kapitole 4.3 vyplývá, že model při zvětšení míry generalizace, při použití 75 % zaměřených bodů, neztrácí na své přesnosti. Ba naopak, v případě kritéria RMSE se model M1_75 jeví jako přesnější. To je však nejspíš zapříčiněno malým počtem verifikačních bodů. Při další generalizaci (na 50 % použitých bodů) už je ztráta přesnosti patrná, avšak nikoliv výrazná. Je třeba také brát v úvahu, že při výběru bodů, které byly odstraněny z původní kompletní sady, nebyl brán ohled na tvar terénu, na výraznější lomové hrany. Tím také mohlo dojít k dalšímu snížení pravdivosti modelu. Model M1_25 s průměrnou absolutní odchylkou 0,8 m se jeví ve srovnání s ostatními jako naprosto nepoužitelný. Vzhledem k algoritmu vyhledávání rozvodnice (kapitola 3.5.5), který postupuje od uzávěrového profilu vždy po spádnicích vzhůru, se mi jako nejvhodnější jeví výrazné zhuštění bodů v blízkosti uzávěrového profilu a zvýšení míry generalizace v horní části povodí. Generalizace terénu v blízkém okolí uzávěrového profilu, způsobuje nepřesnost v celém průběhu rozvodnice.
Okrajové chyby Při triangulaci dochází v okrajových částech modelu k tvorbě trojúhelníků podlouhlého tvaru, což pro reprezentaci terénu není ideální. V těchto trojúhelnících se pak může sklon a orientace svahu značně odlišovat od skutečného terénu. Proto je vhodné zvolit dostatečný přesah měřených dat přes okraje povodí, aby nedošlo k ovlivnění výsledku při hledání průběhu rozvodnice. Dalším problémem při okraji modelu je vytváření hrází, které způsobují tvorbu bezodtokých oblastí značného rozsahu. Tyto oblasti nelze odstranit algoritmem pro běžné
- 62 -
odstraňování bezodtokých oblastí bez ztráty informační hodnoty modelu, proto je potřeba v některých místech modelu vkládat tzv. technické body, které upraví triangulační síť tak, aby byly hráze odstraněny a mohlo docházet k volnému odtoku vody z modelu. Následně pak lze aplikovat algoritmus pro odstranění ostatních bezodtokých oblastí.
Srovnání TIN modelů s modelem na bázi čtvercové sítě Modely TIN a model na bázi čtvercového gridu byly porovnány pouze graficky. Pro model na bázi čtvercové sítě byl použit nejjednodušší a také nejpoužívanější algoritmus pro trasování odtoku a nalezení rozvodnice. Tento algoritmus určuje pro každý čtverec jeden z osmi směrů odtoku, z těchto vypočtených informací, posléze vyhodnocuje tvar povodí. Oproti tomu vytvořený model na bázi nepravidelné trojúhelníkové sítě dovoluje určit pro každý trojúhelník libovolný směr odtoku. Výsledky, které model na bázi TIN dává, jsou proto plynulejší, bez nepatřičných ostrých lomů rozvodnice. Pro srovnávací model na bázi čtvercové sítě by bylo vhodnější použít některý jiný algoritmus pro trasování odtoku, který by více odpovídal podmínkám modelu TIN.
5.3 Analýza povodí Pro každý z vytvořených modelů byla provedena analýza povodí se stejným uzávěrovým bodem. Byla nalezena rozvodnice a zjištěny některé terénní charakteristiky povodí.
Algoritmus pro hledání rozvodnice Algoritmus pro hledání rozvodnice (podrobněji popsaný v kapitole 3.5.5) pracuje na základě hledání průsečíků spádnic s hranami trojúhelníků od uzávěrového profilu na levou i pravou stranu povodí vzhůru až do bodu, který je lokálním maximem modelu. Poté jsou obě větve rozvodnice spojeny do jediného řetězce párů souřadnic, který definuje polygon povodí. V případě povodí Modrava 1 je tento algoritmus plně dostačující, neboť rozvodnice prochází pouze jediným vrcholem – lokálním maximem. V případě povodí s jinou morfologií – s více lokálními maximy – by došlo k hrubé chybě lokalizace rozvodnice. Dvě lokální maxima nejblíže k uzávěrovému profilu by byla spojena přímkou,
- 63 -
body rozvodnice mezi nimi by byly vynechány, čímž by došlo zásadní chybě ve výpočtu plochy povodí a všech ostatních charakteristik. Pro povodí s více vrcholy by musel být algoritmus doplněn o nalezení zbývajících částí rozvodnice. To by mohlo probíhat buďto od zmíněných nejbližších lokálních maxim algoritmem obdobným jako jsou kroky 4 a 5 základního algoritmu (kap. 3.5.5) ovšem opakovaně a s hledáním extrémního sklonu ve směru dolů i nahoru. Další eventualitou je provést hledání sedla mezi lokálními maximy a další část rozvodnice hledat obdobným způsobem jako dosud. Nalezení sedla by však v určitých případech mohlo působit značné problémy.
Charakteristiky povodí V každém z vytvořených modelů byly zjišťovány některé charakteristiky povodí. Ty byly porovnány s doposud udávanými hodnotami, které byly odhadnuty pomocí map. V základním modelu M1 a v modelu M1_50 se plocha povodí po zaokrouhlení s odhadnutou plochou shoduje. Ve zbylých modelech se vypočtená hodnota od odhadované liší v řádech tisíců až desetitisíců metrů čtverečních.(V případě modelu M1_75 je rozdíl činí 5,9 %, v případě modelu M1_25 19,0 %.) V případě průměrného sklonu svahů se všechny čtyři modely shodují v řádu stupňů. Vypočtené hodnoty se nacházejí v rozmezí 6,48 až 6,82 stupně. (Odchylka největšího a nejmenšího modelovaného průměrného skonu svahů činí tedy pouze 0,34°.) Odhadovaný průměrný sklon je 0,09, což odpovídá 5,14°. Rozdíl mezi odhadovanou hodnotou a vypočtenými hodnotami je patrně užitím zjednodušeného vztahu pro odhad charakteristiky.
- 64 -
6 Závěr Cílem předkládané diplomové práce bylo vytvořit digitální model terénu šumavského experimentálního povodí Modrava 1 na bázi nepravidelné trojúhelníkové sítě pomocí vlastního souboru algoritmů. Následně měla být pomocí modelu nalezena orografická rozvodnice a vyhodnoceny základní terénní charakteristiky povodí. Všechny vymezené cíle se v práci podařilo splnit. Výstupem práce je především soubor algoritmů spustitelných v prostředí R sloužící pro vytváření, zobrazování a analýzu digitálního modelu terénu – včetně návodu na jeho obsluhu (kapitola 3.6) – a dále vlastní digitální model povodí Modrava 1 ve čtyřech provedeních podle použité míry generalizace terénu. Všechny vytvořené modely byly vzájemně porovnány. Výstupy analýzy modelů jsou uvedeny v kapitole 4 a doplněny dalšími grafickými přílohami na konci práce. Vytvořený soubor algoritmů je zcela dostačující pro analýzu digitálního modelu povodí Modrava 1, avšak – jak je popsáno v diskusi v kapitole 5.3 – není zcela univerzální. Při jeho použití na povodí s odlišnou morfologií by vedlo k výraznému zkreslení modelované rozvodnice a tím i chybným analýzám všech terénních charakteristik povodí. Dále při ladění programu docházelo k občasným kolapsům vlivem chyb v algoritmizaci některých jednotlivých podfunkcí, všechny odhalené chyby byly opraveny, avšak nelze vyloučit objevení dalších. Vývoj a ladění programu je velmi dlouhodobá záležitost. V budoucnu bych rád program dopracoval do verze, která by byla univerzální – funkční pro libovolnou sadu tachymetrických dat, nedocházelo by v ní ke zmíněné chybě v lokalizaci rozvodnice, a ve které by byly přidané i další aplikace (jako například zobrazení modelu ve formě vrstevnicové mapy, nalezení údolnice povodí a další). Další uvažovanou možností je program přepracovat jako samostatnou aplikaci s vlastním, přístupnějším uživatelským rozhraním. Dopracovaná verze by později mohla sloužit jako názorná pomůcka při výuce hydrologie, případně by mohla posloužit jako základ pro srážko-odtokové modely.
- 65 -
7 Seznam informačních zdrojů BARTÁK, V. (2008): Algoritmy pro zpracování digitálních modelů terénu s aplikacemi v hydrologickém modelování. Diplomová práce, nepublikováno. Česká zemědělská univerzita v Praze, Fakulta životního prostředí, Katedra vodního hospodářství, Praha. BAŠTA, P. (2008): Digitální model terénu povodí Modrava 2. Diplomová práce, nepublikováno. Česká zemědělská univerzita v Praze, Fakulta Životního prostředí, Katedra vodního hospodářství, Praha. BAYER, T. (2013): Rovinné triangulace a jejich využití. Katedra aplikované geoinformatiky
a
kartografie,
Přírodovědecká
fakulta
UK.
Praha.
Online:
http://web.natur.cuni.cz/~bayertom/Adk/adk5.pdf . cit: 2013-04-20. BEYER, W.H. (1987): CRC Standard Mathematical Tables, 28th ed. Boca Raton, FL: CRC Press, str. 123-124. DICKERSON, M.T. – DRYSDALE, R.L.S. – MCELFRESH S.A. – WELZL E. (1997): Fast Greedy Triangulation Algorithms. Computational Geometry, 8 (1998), 67-86. FELDMAN, A.D. – DEVANTIER, B.A. (1993): Review of GIS Applications in Hydrologic Modelling. Journal of Water Resources Planning and Management, Vol. 119., No.2, March/Apríl 1993. Hydrologic Engeneering Center, US Army Corps of Engeneers. GILBERT, P.D. (1979): New Results in Planar Triangulations. Master’s Thesis, University of Illinois, Urbana, 1979. GUDMUNDSSON, J. – HAVERKOT J.H. – VAN KREVELD, M. (2005): Constrained Higher Order Delaunay triangulations. Computational Geometry, vol. 30, Issue 3, March 2005, Pages 271–277. HENGL, T. – GRUBER, S. – SHRESTHA, D.P. (2003): Digital Terrain Analysis in ILWIS. Lecture Notes and User Guide.
- 66 -
Online: https://www.itc.nl/library/Papers_2003/misca/hengl_digital.pdf, cit: 2008-04-01. HJELLE, Ø. – DÆHLEN, M. (2010): Triangulations and Applications. Berlin, Springer. 234 s. HRÁDEK, F. – KUŘÍK, P. (2002): Hydrologie. Česká zemědělská univerzita, Lesnická fakulta. Praha. JONES, J. A. A. – WRIGHT, S. J. – MAIDMENT, D. R. (1990): Watershed Delineation with Triangle-based Terrain Models. Journal of Hydraulic Engeneering 116:1232-51. KLIMÁNEK, M. (2006): Digitální modely terénu. Brno: Mendelova zemědělská a lesnická univerzita v Brně, 2006. 85 s. KOLINGEROVÁ, I. (1999): Rovinné triangulace. Habilitační práce, Západočeská univerzita, Plzeň. KOLINGEROVÁ, I. (2003): On Triangulations. University of West Bohemia, Plzeň. KVHEM (2013): Experimentální povodí Modrava. Katedra vodního hospodářství FŽP ČZU, Praha, online: http://www.kvhem.cz/vyzkum/povodi-modrava/, cit.: 20.3. 2013. LEE, J. (1991): Comparison of Existing Methods for Building Triangullar Irregular Network Models of Terrain from Grid Digital Elevation Models. Int. J. Geographical Information Systems, 1991, Vol. 5, No. 3, 267-285. LEIFER, F. (2006): Delaunayho triangulace a její aplikace. Diplomová práce, nepublikováno. VUT v Brně, FSI ústav automatizace a informatiky. LLOYD, E.L. (1977): On Triangulations of a Set of Points in the Plane. in: Proc. 18th FOCS 228-240. MOORE, I.D. – GRAYSON, R.B. – LADSON, A.R. (1991): Digital Terrain Modeling: A Review of Hydrological, Geomorphological, and Biological Applications. Hydrological Processes, Vol. 5(1991), 3-30. MOORE, I.D. – O‘LOUGHLIN, E.M. – BURCH, G.J. (1988): A Contour-based
- 67 -
Topographic Model for Hydrological and Ecological Applications. Earth Surface Processes and Landforms, vol. 13, 305-320. ONSTAD, C. A. – BRAKENSIEK, D. L. (1968): Watershed Simulation by Stream Path Analogy. Water Resources Research 4, 965-971. PAVLÁSEK, J. – MÁCA, P. – ŘEDINOVÁ J. (2006): Analýza hydrologických dat z modravských povodí. Journal of Hydrology and Hydromechanics, 54 (2006/2), 207-216, Institute of Hydrology of the Slovak Academy of Science. PAVLÁSEK, J. (2013): Ústní sdělení. (2013-04-10). PEUCKER, T.K. – FLOWER, R.J – LITTLE J.J. (1978): The Trainagulated Irregular Network. Department of Geography, Simon Fraser University, Burnaby, B.C., Canada. RENKA, R.J. (1996): Algorithm 751: TRIPACK: a Constrained Two-dimensional Delaunay Triangulation Package. ACM Transactions on Mathematical Software. 22, 1-8. SHOJAEE, D. – HELALI, H. – ALESHEIKH, A. A. (2006): Triangulation for Surface Modelling. In: Ninth International Symposium on the 3D Analysis of Human Movement”, France. SU, P. – DRYSDALE, R.L.S. (1997): A Comparison of Sequential Delaunay Triangulation Algorithms. Computational Geometry, 7 (1997), 361-385. TAJCHMAN, S. J. (1981): On Computing Topographic Characteristic of Mountain Catchment. Canadian Journal of Forest Research 11:768-74. VIVONI E. R. et al. (2004): Generation of Triangulated Irregular Networks Based on Hydrological Similarity. Journal of Hydrologic Engeneering, Vol. 9, No. 4, July 1, 2004. str. 288-302. VOJENSKÝ KARTOGRAFICKÝ ÚSTAV (1991): Šumava – Povydří: turistická mapa 1:50 000. Edice Klubu českých turistů č. 65., 1. vyd., Praha. 548 mm x 742 mm. VÚGTK (2013): Terminologický slovník zeměměřictví a katastru nemovitostí. Výzkumný
- 68 -
ústav geodetický, topografický a kartografický, Terminologická komise ČUZK, Praha. Oniline: http://www.vugtk.cz/slovnik, cit: 2013-04-19. WILSON, J.P. – GALLANT, J.C. (2000): Terrain analysis: Principles and Applications. John Wiley & Sons, 2000. 479 s. YU, S. - VAN KREVELD, M. – SNOEYING J. (1997): Drainage Queries in TINs: from Local to Global and Back Again. In M. J. Kraak and M. Molenaar (eds.), Advances in GIS research II: Proceedings of the Seventh International Symposium on Spatial Data Handling. Lonon, Tailor & Francis, 829-42. ZÁBRANSKÝ, J. (2005): Triangulace povrchů a úlohy na nich. Diplomová práce, nepublikováno. Západočeská univerzita, Fakulta aplikovaných věd, Plzeň.
- 69 -
8 Seznam příloh
Příloha 1
Počet stran
Předmět přílohy
1
Funkce
programu
ArcGIS
použité
pro
tvorbu
srovnávacího modelu na bázi čtvercové sítě Příloha 2
4
Další výstupy modelu M1
Příloha 3
2
Grafické výstupy modelu M1_75
Příloha 4
2
Grafické výstupy modelu M1_50
Příloha 5
2
Grafické výstupy modelu M1_25
Příloha 6
1
Falešná hráz - problematická oblast modelu
Příloha 7
3
Výstupy srovnávacího modelu na bázi čtvercové sítě
Příloha 8
2
Seznam funkcí vytvořeného programu
Příloha 9
2
Obsah CD
Příloha 10
-
CD
- 70 -
Příloha 1, Str. 1/1
Příloha 1: Funkce programu ArcGIS použité pro tvorbu srovnávacího modelu na bázi čtvercové sítě
1.
Natural neighbor
2.
Fill
3.
Flow directions
4.
Slope
5.
Flow accumulation
6.
Watershad
7.
Raster to polygon
8.
Feature to line
Příloha 2, Str.1/4
Příloha 2: Další výstupy modelu M1
Analýza hran modelu M1
Příloha 2, Str.2/4
Tabulky rozložení orientace a sklonu v celém modelu M1 Orientace [°] 0-10 10-20 20-30 30-40 40-50 50-60 60-70 70-80 80-90 90-100 100-110 110-120 120-130 130-140 140-150 150-160 160-170 170-180 180-190 190-200 200-210 210-220 220-230 230-240 240-250 250-260 260-270 270-280 280-290 290-300 300-310 310-320 320-330 330-340 340-350 350-360
Plocha [m2] 2826.03 18380.45 23351.99 25376.02 15097.92 32363.14 58148.72 75187.73 44010.09 19009.23 9479.99 11425.72 9863.31 4660.64 2264.53 49.66 28.60 3391.82 1116.73 0.00 7729.70 2945.94 1731.96 0.00 0.00 631.47 401.74 1812.52 1050.65 1932.31 857.53 1536.01 2538.42 2343.47 2275.71 1912.72
Sklon [°] 0-2 2-4 4-6 6-8 8-10 10-12 12-14 14-16 16-18 18-20 20-22 22-24 24-26 26-28 28-30 30-32 32-34 34-36 36-38 38-40 >40
Plocha [m2] 14271.55 57082.47 141360.30 104964.00 46398.05 10174.35 6984.90 1186.97 232.70 648.46 720.65 295.53 111.82 759.18 0.58 25.77 44.87 0.00 29.92 0.00 440.43
Příloha 2, Str.3/4
Ukázka znázornění tras odtoku v modelu funkcí placeDrop. Pro ilustraci byly náhodně vybrány body uvnitř i vně povodí.
Příloha 2, Str.4/4
Ukázky třírozměrného znázornění modelu. Vzhledem k malé členitosti terénu je osa Z v trojnásobném měřítku.
Příloha 3, Str.1/2
Příloha 3: Grafické výstupy modelu M1_75
Grafické znázornění orientace svahů, směrů odtoku a průběhu rozvodnice modelu M1_75
Příloha 3, Str.2/2
Grafické znázornění sklonu svahů a průběhu rozvodnice modelu M1_75
Příloha 4, Str.1/2
Příloha 4: Grafické výstupy modelu M1_50
Grafické znázornění orientace svahů, směrů odtoku a průběhu rozvodnice modelu M1_50
Příloha 4, Str.2/2
Grafické znázornění sklonu svahů a průběhu rozvodnice modelu M1_50
Příloha 5, Str.1/2
Příloha 5: Grafické výstupy modelu M1_25
Grafické znázornění orientace svahů, směrů odtoku a průběhu rozvodnice modelu M1_25
Příloha 5, Str.2/2
Grafické znázornění sklonu svahů a průběhu rozvodnice modelu M1_25
Příloha 6, Str. 1/1
Příloha 6: Falešná hráz - problematická oblast modelu
Falešná hráz – na obrázku znázorněna červeně (červené trojúhelníky - orientované k jihu) na severním okraji modelu. Technické body byly do modelu přidány následně pro odstranění chyby.
Příloha 7, Str. 1/3
Příloha 7: Grafické výstupy srovnávacího modelu na bázi čtvercové sítě
Směry odtoku určené modelem na bázi čtvercové sítě pomocí programu ArcGIS 10
Příloha 7, Str. 2/3
Sklony svahu určené modelem na bázi čtvercové sítě pomocí programu ArcGIS 10
Příloha 7, Str. 3/3
Srovnání rozvodnic vytvořených TIN modelem M1 a modelem na bázi gridu na pozadí modelu TIN. Jsou patrné značné odlišnosti ve tvaru povodí.
Příloha 8, Str.1/2
Příloha 8: Seznam funkcí vytvořeného programu basinAnalysis
Obalová funkce, provede analýzu povodí se zadaným uzávěrovým profilem (nalezení rozvodnice, plocha povodí, střední sklon svahů, celková orientace povodí).
basinDivide
Nalezne rozvodnici povodí se zadaným uzávěrovým profilem.
componentArea
Počítá plochu, kterou trojúhelník přispívá do povodí.
createArcsTable
Vytváří seznam hran a hrany analyzuje.
createMesh
Vytváří triangulaci nad zadanou množinou bodů.
distance
Interaktivní funkce, měří vzdálenost dvou bodů v modelu.
drainlessRemove
Vyhledává a odstraňuje bezodtoké oblasti.
findNormalVector
Najde normálové vektory všech trojúhelníků.
findTDirection
Vypočte orientaci všech trojúhelníků.
findTSlope
Vypočte sklon všech trojúhelníků.
findZ
Lineárně interpoluje souřadnici Z libovolného bodu ze souřadnic vrcholů trojúhelníku v němž leží.
flow
Najde trasu odtoku ze zadaného bodu až k okraji modelu.
flowaross
Hledá další bod trasy odtoku v případě tečení přes trojúhelník.
flowalong
Hledá další bod trasy odtoku v případě tečení údolím podel trojúhelníku.
ftIntersections
Nachází průsečíky spádnice trojúhelníku procházející zadaným bodem s hranami trojúhelníku ve kterém leží.
globalize
Globalizuje lokální proměnné při tvorbě modelu.
intersection
Nachází průsečík dvou přímek zadaných bodem a směrovým vektorem.
loaddata
Načítá zdrojová data z textového souboru.
mycolorpalette
Vytváří barevnou paletu pro zobrazování orientace svahů.
mycolorpalette2
Vytváří barevnou paletu pro kombinované zobrazování orientace a sklonu svahů.
nearestNode
Interaktivní funkce. Nachází nejbližší uzlový bod a informace o něm vypisuje do konzole.
Příloha 8, Str.2/2
placeDrop
Interaktivní funkce. Zobrazuje trasu odtoku ze zadaného místa.
plDistance
Počítá vzdálenost bodu od přímky zadané bodem a směrovým vektorem.
plotBasinDivide
Zobrazuje rozvodnici do aktivního grafického okna.
plotModel2D
Vytváří dvourozměrné grafické zobrazení modelu.
plotModel3D
Vytváří třírozměrné grafické zobrazení modelu v okně rgl.
plotRidgeEtThalweg Zobrazuje lokální hřbetnice a údolnice do modelu v aktivním okně. saveP
Ukládá projekt.
sdDistribution
Vypočítává rozložení sklonu nebo svahů v povodí.
siftacc
Vytřídí nepoužitelné body z hlediska přesnosti udávané měřicím přístrojem.
siftdist
Vytřídí zdvojené body z hlediska jejich vzájemné vzdáenosti.
sourceAll
Načte všechny potřebné funkce do paměti programu R.
thisPoint
Interaktivní funkce. Do konzole zobrazí informace o libovolném bodu na povrchu modelu.
TINmodel
Obalová funkce. Spouští sekvenci dalších funkcí, jejichž výsledkem je model připravený k analýze.
triangleArea
Počítá plochu trojúhelníku.
verification
Provádí verifikaci modelu z hlediska kritérií MAE a RMSE.
whichTriangle
Zjišťuje, kterému trojúhelníku náleží bod o zadaných souřadnicích.
zoomIN
Interaktivní funkce. Zobrazuje výřez modelu v novém grafickém okně.
Příloha 9, Str.1/2
Příloha 9: Obsah CD DP_Herout.pdf
Text diplomové práce
WorkingDirctory
Pracovní adresář basinAnalysis.R
Soubory R - zdrojové kódy algoritmů
basinDivide.R
…
componentArea.R
…
createArcsTable.R
…
createMesh.R
…
distance.R
…
drainlessRemove.R
…
findNormalVector.R
…
findTDirection.R
…
findTSlope.R
…
findZ.R
…
flow.R
…
flowaross.R
…
flowalong.R
…
ftIntersections.R
…
intersection.R
…
loaddata.R
…
mycolorpalette.R
…
mycolorpalette2.R
…
nearestNode.R
…
placeDrop.R
…
plDistance.R
…
plotBasinDivide.R
…
plotModel2D.R
…
plotModel3D.R
…
plotRidgeEtThalweg.R
…
saveP.R
…
sdDistribution.R
…
Příloha 9, Str.2/2
siftacc.R
…
siftdist.R
…
sourceAll.R
…
thisPoint.R
…
TINmodel.R
…
triangleArea.R
…
verification.R
…
whichTriangle.R
…
zoomIN.R
…
ProjectM1
Uložený projekt (M1)
ProjectM1_25
… (Model M1_25)
ProjectM1_50
… (Model M1_50)
ProjectM1_75
… (Model M1_75)
Data.txt
Textový soubor – zdrojová data (M1)
Data_25.txt
… (Model M1_25)
Data_50.txt
… (Model M1_50)
Data_75.txt
… (Model M1_75)
Data_bezTB.txt
… (původní - bez technických bodů)
VerifikaceM1.txt
Textový soubor – výsledky validace (M1)
VerifikaceM1_25.txt
… (Model M1_25)
VerifikaceM1_50.txt
… (Model M1_50)
VerifikaceM1_75.txt
… (Model M1_75)
Packages
Adresář s externími balíčky