Metoda SCS CN

Teoretický popis metody, základní symboly použité v dalším textu, popis vstupních dat a navrhovaný postup je součástí školení GRASS GIS pro pokročilé.

Postup zpracovaní v QGIS

Krok 1

Sjednocení hlavních půdních jednotek a komplexního průzkumu půd

V prvním kroku založíme projekt a pomocí ^{Přidaní vektorové vrstvy} a ^{Přidat vrstvu s odděleným textem} vložíme do panelu vrstev vstupní data - vektorová data ve formátu Esri Shapefile hpj.shp, kpp.shp, landuse.shp, povodi.shp a soubory s odděleným textem hpj_hydrsk.csv, kpp_hydrsk.csv, lu_hydrsk_cn.csv, viz školení QGIS pro začátečníky.

Následně sjednotíme vrstvu hlavních půdních jednotek (hpj) a komplexního průzkumu půd (kpp). Využijeme nástroj geoprocessingu ^Sjednotit, který najdeme v záložce Vektor ‣ Nástroje geoprocessingu. Vznikne nová vektorová vrstva hpj_kpp.

Krok 2

Připojení informací o hydrologické skupině

Tabulku hpj_hydrsk můžeme připojit k atributům nové vektorové vrstvy hpj_kpp pomocí společného sloupce HPJ (obr. 1). Pravým tlačítkem myši v panelu vrstev u hpj_kpp zvolíme Vlastnosti a v dialogovém okně přejdeme do záložky ^Připojení. Kliknutím na ikonku spustíme dialogové okno s nastavením pro připojení (obr. 2), viz kapitola Připojení tabulkových dat ve školení QGIS pro začátečníky.

Obrázek 1: Společný atribut HPJ a hydrologické skupiny hlavních půdních jednotek.

Obrázek 2: Připojení tabulky k vektorové vrstvě na zákládě společného atributu.

Tímto způsobem připojíme tabulky s informacemi o hydrologických skupinách (obr. 3).

Obrázek 3: Zobrazení připojených vrstev ve vlastnostech vrstvy hpj_kpp.

Poznámka

V některých případech nemusí připojení tabulek proběhnout korektně, např. pokud se liší datový typ sloupečků použitých pro připojení. Typicky textové vs. číselné pole, např. '09' vs. '9'.

Poté otevřeme atributovou tabulku hpj_kpp, zapneme editační mód ikonkou ^{Přepnout editaci} a pomocí kalkulačky polí ^{Otevřít kalkulačku polí} vytvoříme nový atribut. Použijeme připojené atributy o hydrologické skupině (hpj_HydrSk z hlavních půdních jednotek a kpp_HydrSk z komplexního průzkumu půd). Primárně použijeme hydrologickou skupinu pro hlavní půdní jednotky. Kde informace není dosupná - hodnota NULL, tam použijeme kpp_HydrSk (obr. 4) a výsledek znázorníme (obr. 5).

CASE WHEN "hpj_HydrSk" IS NULL THEN "kpp_HydrSk" ELSE "hpj_HydrSk" END

Obrázek 4: Vytvoření atributu s informacemi o hydrologické skupině pro elementární plochy.

Obrázek 5: Hydrologické skupiny elementárních ploch v zájmovém území.

Při pohledu na legendu na obr. 5 je možno si všimnout, že kódy hydrologických skupin jako (A)B, A(B), AB a podobně by bylo vhodné sjednotit. K tomu použijeme editační mód a atributové dotazy. V hlavní liště anebo v liště atributové tabulky zvolíme volbu ^{Vybrat prvky pomocí vzorce} pomocí které vybereme elementární plochy s hydrologickou skupinou (A)B a A(B), potom zapneme editační režim, spustíme ^{Kalkulačka polí} a aktualizujeme existujíce atributy hydrsk vybraných prvků (obr. 6). Obdobně postupujeme pro další kódy. Výsledek je prezentován na obr. 7.

Obrázek 6: Sjednocení hodnot atributů pomocí kalkulátoru polí.

Poznámka

Pro sjednocení hodnot je možno použít také výraz CASE:

CASE WHEN "hydrsk"  =  'B(C)' THEN replace("hydrsk",'B(C)','BC') ELSE "hydrsk" END

a

CASE WHEN "hydrsk"  =  'C(D)' THEN replace("hydrsk",'C(D)','CD') ELSE "hydrsk" END

Obrázek 7: Sjednocené hydrologické skupiny elementárních ploch v zájmovém území.

Do této fáze je možné používat QGIS relativně bez problémů. Dále však budeme přidávat informace o využití území pro každou elementární plochu pomocí operace průniku. Při větších objemech dat mohou být nástroje geoprocessingu časově náročné a nestabilní. Pro další řešení tedy použijeme výpočetně stabilnější nástroje systému GRASS GIS.

Tip

Více o systému GRASS v rámci školení GRASS GIS pro začátečníky.

Vytvoření lokace a mapsetu

Data, ke kterým GRASS přistupuje, udržuje v pevné 3-úrovňové struktuře (databáze, lokace a mapset), viz Struktura dat - koncept lokací a mapsetů ze školení GRASS GIS pro začátečníky. Z hlavní lišty menu vybereme Zásuvné moduly ‣ GRASS ‣ Nový mapset.

Obrázek 8: Zásuvný modul GRASS - vytvoření nového mapsetu.

V dialogovém okně se objeví předvolená cesta k hlavnímu adresáři databáze GRASS (obvykle adresář s názvem grassdata). V případě, že tento adresář obsahuje již nějaké lokace, vybereme tu, ve které chceme pracovat anebo si vytvoříme novou. Nastavíme souřadnicový systém a výpočetní region (viz školení GRASS GIS pro začátečníky). Kromě mapsetu PERMANENT, který se vytvoří automaticky, je vhodné zadat i název nového mapsetu, ve kterém budou probíhat výpočty. Mapset se automaticky otevře jako aktuální mapset. V záložce Region dialogového okna nástrojů GRASS je možné měnit rozsah výpočetní oblasti výběrem v mapovém okně QGIS pomocí Select the extent by dragging on canvas (obr. 9). Zároveň je zde možno nastavit prostorové rozlišení.

Obrázek 9: Vytvoření lokace a mapsetu, nastavení výpočetní oblasti a prostorového rozlišení.

Krok 3

Průnik vrstvy hydrologických skupin s vrstvou využití území

Zájmové území potřebujeme rozdělit na více elementárních ploch. Vrstvy hpj_kpp a landuse, pro které vytvoříme průnik, musíme nejprve naimportovat do mapsetu. Import dat zajišťuje více nástrojů, tzv. modulů (obr. 10). Použijeme například modul v.in.ogr.qgis, který umožňuje načítat vrstvy (jakoby) z prostředí QGIS. Názvy vrstev zachováme stejné.

Obrázek 10: Možnosti importu vektorových vrstev do GRASS mapsetu v prostředí QGIS.

Pokud chceme oveřit, zda se zadané vrstvy po importu v mapsetu nacházejí, použijeme shell. Kliknutím na ^{GRASS shell} spustíme příkazový řádek. Obsah konkrétního mapsetu vypiše modul g.list. Pro výpis vektorových vrstev v aktuálním mapsetu zadáme g.list type=vector mapset=.. Pokud zadáme pouze příkaz g.list, tak se otevře dialogové okno modulu a parametry můžeme zadat interaktivně.

Poznámka

Dokumentaci a povinné parametry každého modulu lze zobrazit zadáním příkazu man před název modulu, například man g.list.

Operaci překrytí, resp. určení průniku vektorových vrstev, zajišťuje modul ^{v.overlay.and}, který spustíme z Vektor ‣ Prostorová analýza ‣ Překrytí (obr. 11). Výslednou vrstvu průniku nazveme hpj_kpp_landuse.

Poznámka

Počet záznamů v atributové tabulce se průnikem prvků výrazně zvýší. Což je zapříčiněno hlavně tím, že QGIS zobrazuje záznamy pro multiprvky jako duplicitní.

Obrázek 11: Modul pro určení průniku dvou vektorových vrstev.

Tip

V příkazovém řádku můžeme vypsat například:

• seznam tabulek v aktuálním mapsetu, resp. jejich názvy: db.tables
• seznam atributů konkrétní tabulky: db.columns table=NAZEVTABULKY
• počet záznamů v tabulce: db.select sql='select count(*) from NAZEVTABULKY'

Příklad použití GRASS shell je znázorněn na obr. 12. Modul ^v.db.select vypíšeme hodnoty atributů, modulem ^{v.db.select.where} je možné zadat i podmínku.

Obrázek 12: Zobrazení tabulek a záznamů v příkazovém řádku.

Modul v.out.ogr umožňuje exportovat atributovou tabulku do různých formátů a dále s ni pracovat. Na obr. 13 je znázorněn export do bežného formátu CSV.

Obrázek 13: Export atributů do formátu CSV.

Krok 4

Připojení hodnot odtokové křivky \(CN\)

V dalším kroku je potřeba vytvořit atribut, který bude obsahovat údaje o využití území a o hydrologické skupině půdy dané elementární plochy ve tvaru VyužitíÚzemí_HydrologickáSkupina.

Vytvoříme nový atribut pomocí modulu ^{v.db.add.column}, který nazveme landuse_hydrsk (obr. 14). Potom doplníme hodnoty atributu s využitím modulu ^{v.db.update_op} jako výsledek operace v rámci jedné atributové tabulky. Hodnotu zadáme ve tvaru b_LandUse||'_'||a_hydrsk.

Obrázek 14: Přidání atributu do atributové tabulky s datovým typem text.

Poznámka

Výsledek můžeme zkontrolovat v příkazovém řádku zadaním

v.db.select map=hpj_kpp_landuse columns=cat,b_LandUse,a_hydrsk,landuse_hydrsk where=cat=1

cat|b_LandUse|a_hydrsk|landuse_hydrsk
1|OP|B|OP_B

Dále do mapsetu modulem db.in.ogr importujeme tabulku s hodnotami CN. Nazveme ji lu_hydrsk_cn.

Následně použijeme modul ^v.db.join pomocí kterého připojíme importovanou tabulku k vektorové vrstvě hpj_kpp_landuse a to kvůli přiřazení hodnot CN ke každé elementární ploše řešeného území, viz. obr. 15. Obsah výsledné tabulky je možno oveřit v příkazovém řádku pomocí v.db.select map=hpj_kpp_landuse where=cat=1.

Důležité

Jednotlivé atributy v tabulkách, které spojujeme , nesmí obsahovat stejné názvy atributů. Tento problém lze vyřešit zavoláním modulu v.db.join z GUI systému GRASS a volbou subset_columns, která v GRASS pluginu chybí.

Obrázek 15: Připojení tabulky k existující atributové tabulce vektorové vrstvy.

Poznámka

Tento způsob spojení atributových dat je alternativou k operaci záložky ^Připojení ve vlastnostech vektorové vrstvy, viz Krok 2.

Krok 5

Sjednocení průniku vrstvy hydrologických skupin a využití území s vrstvou povodí

Hodnoty návrhových srážek s různou dobou opakovaní do vrstvy přidáme pomocí modulu ^v.overlay.or. Sjednocení předchází import vrstvy povodí s informacemi o srážkách do mapsetu, přičemž postup je obdobný jako při importu vektorových vrstev v úvodní části.

Ukázka záznamu (vybrané atributy) atributové tabulky nově vytvořené vektorové vrstvy hpj_kpp_lu_pov pro 2-letý úhrn srážek v mm s dobou trvaní 120 min:

v.db.select map=hpj_kpp_lu_pov columns=cat,a_CN,b_H_002_120 where"cat=1"

cat|a_CN|b_H_002_120
1|80|21.6804582207

Přehled o tom, jak se změnil počet plošných prvků ve vrstvě hpj_kpp_landuse po sjednocení s vrstvou povodí, dostaneme jako výstup modulu v.info, viz. Vektor ‣ Zprávy a statistiky. Standardní zobrazení je uvedeno obr. 16.

Obrázek 16: Výpis základních informací o vektorové vrstvě pomocí modulu v.info.

Tip

Z příkazového řádku je možno spustit nativní grafické uživatelké rozhraní systému GRASS příkazem g.gui. Taktéž je možné zapnout mapové okno (příkaz d.mon), vykreslit v něm konkrétní rastrovou (d.rast) anebo vektorovou (d.vect) vrstvu, přidat měřítko (d.barscale) či legendu (d.legend). Příkazem d.rast.leg vykreslíme rastrovou vrstvu i s legendou.

Dále budeme pracovat především s hodnotami CN. Pro další operace je potřeba, aby typ tohoto atributu byl číselný, na což použijeme funkci cast(). Vytvoříme tedy nový atribut CN s datovým typem integer.

Poznámka pro pokročilé

Vektorovou vrstvu hpj_kpp_landuse je možno převést na rastrovou vrstvu s hodnotami CN a zobrazit v mapovém okně systému GRASS. Začneme vytvořením nového atributu typu integer (modul v.db.addcolumn), pokračujeme jeho editací v.db.update_op a následně spustíme modul ^{v.to.rast.attr}, viz. obr. 17. Pomocí příkazů d.mon wx0, d.rast.leg cn, d.barscale a d.vect povodi type=boundary zobrazíme mapu s CN včetně měřítka legendy v překrytí s vektorovou vrstvou povodí.

Obrázek 17: Konverze vektorové mapy na rastrovou na základě atributu.

Krok 6

Výpočet výměry elementárních ploch, parametru \(A\) a parametru \(I_a\)

Pro každou elementární plochu vypočítame její výměru, parametr \(A\) a \(I_a\).

\[A = 25.4 \times (\frac{1000}{CN} - 10)\]

\[I_a = 0.2 \times A\]

Do atributové tabulky hpj_kpp_lu_pov přidáme nové atributy typu double, konkrétně vymera, A, I_a. Poté vypočítame jejich příslušné hodnoty. Postupujeme obdobně jako při tvorbě atributu s hodnotami o využití území a hydrologické skupině (landuse_hydrsk), přičemž pro jejich výpočet použijeme matematické operáce jako sčítaní, odčítaní, násobení a podobně (obr. 18 a obr. 19). Pro určení plochy každé elementární plochy využijeme modul z kategorie Vektor ‣ Zprávy a statistiky, modul ^v.to.db.

Obrázek 18: Vytvoření více atributů najednou s využitím v.db.addcolumn.

Obrázek 19: Výpočet výměry modulem v.to.db a parametru A modulem v.db.update_op.

Poznámka pro pokročilé

V příkazovém řádku by tyto kroky vypadaly následovně:

v.db.addcolumn map=hpj_kpp_lu_pov columns="vymera double,A double,Ia double"
v.to.db map=hpj_kpp_lu_pov option=area columns=vymera
v.db.update map=hpj_kpp_lu_pov column=A value="24.5 * (1000 / CN - 10)"
v.db.update map=hpj_kpp_lu_pov column=I_a value="0.2 * A"

Krok 7

Výpočet parametru \(H_o\) a parametru \(O_p\) pro každou elementární plochu

Znázornění vektorové vrstvy povodí s návrhovými srážkami v prostředí QGIS je uvedeno na obr. 20 (maximální hodnota atributů H_002_120 představuje 23 mm). Histogramy je možné vykreslit v záložce ^Diagramy ve vlastnostech konkrétní vrstvy.

Obrázek 20: Zobrazení povodí IV. řádu s návrhovými srážkami.

Vypočítáme výšku přímého odtoku v mm jako parametr \(H_o\) a objem jako parametr \(O_{p}\).

\[H_O = \frac{(H_S − 0.2 \times A)^2}{H_S + 0.8 \times A}\]\[O_P = P_P \times \frac{H_O}{1000}\]

V dalších krocích budeme uvažovat průměrný úhrn návrhové srážky \(H_{s}\) = 32 mm. Při úhrnu s dobou opakovaní 2 roky (atribut H_002_120) či dobou 5, 10, 20, 50 anebo 100 roků by byl postup obdobný.

Důležité

Hodnota v čitateli vztahu pro \(H_o\) musí být kladná, resp. nelze umocňovat záporné číslo. V případě, že čitatel je záporný, výška přímého odtoku je rovná nule. Pomůžeme si novým atributem v atributové tabulce, který nazveme HOklad.

Postupujeme obdobně jako na obr. 18 a obr. 19 anebo pomocí příkazového řádku.

v.db.addcolumn map=hpj_kpp_lu_pov columns="HOklad double, HO double, OP double"
v.db.update map=hpj_kpp_lu_pov column=HOklad value="(32 - (0.2 * A))"

Záporným hodnotám HOklad přiřadíme konstantu 0 modulem ^{v.db.update_query} (obr. 21). Atributy HO a OP vyplníme modulem ^{v.db.update_op}.

v.db.update map=hpj_kpp_lu_pov column=HO value='(HOklad * HOklad)/(32 + (0.8 * A))'
v.db.update map=hpj_kpp_lu_pov column=OP value="vymera * (HO / 1000)"

Obrázek 21: Přiřazení konstatní hodnoty atributu v případě splnění podmínkt dotazu modulem v.db.update_query.

Tip

Přiřazení konstanty 0 pro záporné HOklad je možno zkontrolovat tak jako je prezentovano na obr. 22.

Obrázek 22: Kontrola editace záporných hodnot v příkazovém řádku.

Krok 8

Vytvoření rastrových vrstev výšky a objemu přímého odtoku

Modulem ^{v.to.rast.attr} vytvoříme z vektorové vrstvy hpj_kpp_lu_pov rastrové vrstvy ho a op. Výsledky vizualizované v prostředí QGIS jsou uvedeny na obr. 23.

Varování

Před samotnou rasterizací je nutné korektně nastavit výpočetní region.

Obrázek 23: Zobrazení výšky a objemu přímého odtoku pro elementární plochy v prostředí QGIS.

Krok 9

Výpočet průměrných hodnot výšky a objemu přímého odtoku pro povodí

V dalším kroku vypočítáme průměrné hodnoty přímého odtoku pro každé povodí v řešeném území. Modul ^v.rast.stats počítá základní statistické informace rastrové vrstvy na základě vektorové vrstvy a ty ukladá do nových atributů v atributové tabulce. Dialogové okno je uvedeno na obr. 24.

Obrázek 24: Dialogové okno modulu v.rast.stats.

Vektorovou vrstvu povodí potom převedeme do podoby rastrové vrstvy, přičem jako klíčový atribut použime ho_average, resp. op_average. Výstup zobrazený v prostředí QGIS je na obr. 25.

Obrázek 25: Zobrazení průměrné výšky a objemu přímého odtoku pro povodí v prostředí QGIS.