153YZOD Zpracování obrazových dat 2006 - 4. cvičení: Porovnání verzí
m (→rastrová algebra: korekce) |
m (upozorneni) |
||
(Není zobrazeno 15 mezilehlých verzí od 2 dalších uživatelů.) | |||
Řádek 1: | Řádek 1: | ||
<div class="boilerplate metadata" id="attention" style="border: 3px double #c32c2c; margin-bottom: 2em; padding: 0 .25em; background-color: #f96666;"> | |||
<p style="font-size: 150%; text-align: center;"><br/>'''Tento text se vztahuje k akademickému roku 2006/2007. Aktuální stránky cvičení jsou dostupné [[153YZOD Zpracování obrazových dat - cvičení 4|zde]].'''<br/> </p> | |||
<div style="clear: both;"><br style="clear: both; display: none;" /></div> | |||
</div> | |||
{{ZOD}} | |||
[ [[Zpracování obrazových dat]] ] | |||
= Mapová algebra, sčítání, odčítání rastrových dat, podíl obrazu, NDVI = | = Mapová algebra, sčítání, odčítání rastrových dat, podíl obrazu, NDVI = | ||
__TOC__ | __TOC__ | ||
Řádek 42: | Řádek 49: | ||
pravidel pro vlastní tabulku barev v modulu <tt>r.colors</tt>. | pravidel pro vlastní tabulku barev v modulu <tt>r.colors</tt>. | ||
#interaktivní mód modulu r.mapcalc | |||
#interaktivní mód modulu r.mapcalc | # | ||
# | GRASS:~ > r.mapcalc | ||
GRASS:~ > r.mapcalc | |||
Enter expressions, "end" when done. | Enter expressions, "end" when done. | ||
mapcalc> p12 = tm1 / tm2 | mapcalc> p12 = tm1 / tm2 | ||
mapcalc> p23 = tm2 / tm3 | mapcalc> p23 = tm2 / tm3 | ||
mapcalc> end | mapcalc> end | ||
Výhodnější je však předat dané výrazy modulu jako parametr (uzavřít do <u>uvozovek</u>!), nebo | Výhodnější je však předat dané výrazy modulu jako parametr (uzavřít do <u>uvozovek</u>!), nebo | ||
je uložit do textového souboru (například při opětovném použití) | je uložit do textového souboru (například při opětovném použití) | ||
a přesměrovat do modulu přes standardní vstup (podobně jako tabulku barev | a přesměrovat do modulu přes standardní vstup (podobně jako tabulku barev | ||
do <tt>r.colors</tt>). Např. ([http:// | do <tt>r.colors</tt>). Např. ([http://gama.fsv.cvut.cz/~landa/geowikicz_data/zod/cv4/podily.txt podily.txt]): | ||
#zadání rovnice jako parametr | |||
#zadání rovnice jako parametr | # | ||
# | GRASS:~ > r.mapcalc 'p12 = tm1 / tm2' | ||
GRASS:~ > r.mapcalc 'p12 = tm1 / tm2' | GRASS:~ > r.mapcalc 'p23 = tm2 / tm3' | ||
GRASS:~ > r.mapcalc 'p23 = tm2 / tm3' | # | ||
# | #rovnice uložené v textovém souboru | ||
#rovnice uložené v textovém souboru | # | ||
# | GRASS:~ > r.mapcalc < podily.txt | ||
GRASS:~ > r.mapcalc < podily.txt | |||
My bude upřednostňovat způsob zadávání výrazů jako parametru. | My bude upřednostňovat způsob zadávání výrazů jako parametru. | ||
Řádek 75: | Řádek 79: | ||
Obecně výraz pro <tt>r.mapcalc</tt> vypadá následovně: | Obecně výraz pro <tt>r.mapcalc</tt> vypadá následovně: | ||
nova_mapa = vyraz (mapa1, mapa2, ...) | nova_mapa = vyraz (mapa1, mapa2, ...) | ||
kde <tt>vyraz</tt> může zahrnovat názvy rastrových mapových vrstev, celočíselné či reálné | kde <tt>vyraz</tt> může zahrnovat názvy rastrových mapových vrstev, celočíselné či reálné | ||
Řádek 380: | Řádek 382: | ||
naši pozornost si rozhodně zaslouží... | naši pozornost si rozhodně zaslouží... | ||
== podíl | == podíl obrazu == | ||
Podíl obrazu je základní a velmi jednoduchou algebraickou metodou používanou pro detekci příznaků, | Podíl obrazu je základní a velmi jednoduchou algebraickou metodou používanou pro detekci příznaků, | ||
Řádek 391: | Řádek 393: | ||
Jako příklad můžeme vypočítat podíl sedmého a čtvrtého kanálu družicové scény LandSat-TM5: | Jako příklad můžeme vypočítat podíl sedmého a čtvrtého kanálu družicové scény LandSat-TM5: | ||
#jistota... | |||
#jistota... | # | ||
# | GRASS:~ > g.region rast=tm4 | ||
GRASS:~ > g.region rast=tm4 | # | ||
#podíl obrazových kanálů sedm a čtyři (výsledná mapa DCELL) | #podíl obrazových kanálů sedm a čtyři (výsledná mapa DCELL) | ||
# | # | ||
GRASS:~ > r.mapcalc 'p74 = 1.0 * tm7 / tm4' | GRASS:~ > r.mapcalc 'p74 = 1.0 * tm7 / tm4' | ||
GRASS:~ > d.rast p74 | GRASS:~ > d.rast p74 | ||
[[Soubor:ZOD-cv4-p74.png|frame|center|Obr č.3: Podíl obrazu (tm7/tm4) - tabulka barev rainbow]] | [[Soubor:ZOD-cv4-p74.png|frame|center|Obr č.3: Podíl obrazu (tm7/tm4) - tabulka barev rainbow]] | ||
Řádek 406: | Řádek 407: | ||
Příkaz pro tvorbu mapy s NDVI vypadá následovně: | Příkaz pro tvorbu mapy s NDVI vypadá následovně: | ||
#výpočet NDVI | |||
#výpočet NDVI | # | ||
# | GRASS:~ > r.mapcalc 'ndvi = 1.0 * (tm4 - tm3) / (tm4 + tm3)' | ||
GRASS:~ > r.mapcalc 'ndvi = 1.0 * (tm4 - tm3) / (tm4 + tm3)' | |||
Dominuje-li DH v blízkém infračerveném oboru (<u>tm4</u>) nad hodnotou | Dominuje-li DH v blízkém infračerveném oboru (<u>tm4</u>) nad hodnotou | ||
Řádek 422: | Řádek 421: | ||
Pro lepší vizuální interpretaci obrazu zvolíme vlastní tabulku barev | Pro lepší vizuální interpretaci obrazu zvolíme vlastní tabulku barev | ||
([http:// | ([http://gama.fsv.cvut.cz/~landa/geowikicz_data/zod/cv4/tb-ndvi.txt tb-ndvi.txt]), nejprve zjistíme rozsah dat | ||
(<tt>r.info</tt>), pro názornější představu můžeme vytisknout | (<tt>r.info</tt>), pro názornější představu můžeme vytisknout | ||
statistiku dat (<tt>r.report</tt>): | statistiku dat (<tt>r.report</tt>): | ||
#rozsah hodnot | |||
#rozsah hodnot | # | ||
# | GRASS:~ > r.info ndvi | ||
GRASS:~ > r.info ndvi | |||
+----------------------------------------------------------------------------+ | |||
+----------------------------------------------------------------------------+ | | Layer: ndvi Date: Sat Oct 29 14:42:15 2005 | | ||
| Layer: ndvi Date: Sat Oct 29 14:42:15 2005 | | | Mapset: student Login of Creator: martin | | ||
| Mapset: student Login of Creator: martin | | | Location: sevcech | | ||
| Location: sevcech | | | DataBase: /home/martin/grassdata | | ||
| DataBase: /home/martin/grassdata | | | Title: ( ndvi ) | | ||
| Title: ( ndvi ) | | |----------------------------------------------------------------------------| | ||
|----------------------------------------------------------------------------| | | | | ||
| | | | Type of Map: cell Number of Categories: 255 | | ||
| Type of Map: cell Number of Categories: 255 | | | Data Type: DCELL | | ||
| Data Type: DCELL | | | Rows: 1661 | | ||
| Rows: 1661 | | | Columns: 2223 | | ||
| Columns: 2223 | | | Total Cells: 3692403 | | ||
| Total Cells: 3692403 | | | Projection: krovak (zone 0) | | ||
| Projection: krovak (zone 0) | | | N: -957500 S: -1007318 Res: 29.99277544 | | ||
| N: -957500 S: -1007318 Res: 29.99277544 | | | E: -763855.0602047 W: -830529 Res: 29.99277544 | | ||
| E: -763855.0602047 W: -830529 Res: 29.99277544 | | | Range of data: min = -1.000000 max = 0.861314 | | ||
| Range of data: min = -1.000000 max = 0.861314 | | | | | ||
| | | | Data Source: | | ||
| Data Source: | | | | | ||
| | | | | | ||
| | | | | | ||
| | | | Data Description: | | ||
| Data Description: | | | generated by r.mapcalc | | ||
| generated by r.mapcalc | | | | | ||
| | | | Comments: | | ||
| Comments: | | | 1.000000 * (tm4 - tm3) / (tm4 + tm3) | | ||
| 1.000000 * (tm4 - tm3) / (tm4 + tm3) | | | | | ||
| | | +----------------------------------------------------------------------------+ | ||
+----------------------------------------------------------------------------+ | |||
Rastrové hodnoty se tedy pohybují v intervalu od [-1.0; +0.87]. Všimněme si navíc | Rastrové hodnoty se tedy pohybují v intervalu od [-1.0; +0.87]. Všimněme si navíc | ||
typu dat - DCELL (rastrová data s plovoucí desetinnou čárkou). V komentáři (Comments) | typu dat - DCELL (rastrová data s plovoucí desetinnou čárkou). V komentáři (''Comments'') | ||
je navíc uveden výraz, na základě něhož tato data vznikla. | je navíc uveden výraz, na základě něhož tato data vznikla. | ||
#vytisknout statistiku (počet buněk, procenta, plocha v hektarech); přerozdělit na 25 intervalů | |||
#vytisknout statistiku (počet buněk, procenta, plocha v hektarech); přerozdělit na 25 intervalů | # | ||
# | GRASS:~ > r.report -n map=ndvi units=c,p,h nstep=25 | ||
GRASS:~ > r.report -n map=ndvi units=c,p,h nstep=25 | |||
+-----------------------------------------------------------------------------+ | |||
+-----------------------------------------------------------------------------+ | | RASTER MAP CATEGORY REPORT | | ||
| RASTER MAP CATEGORY REPORT | | |LOCATION: sevcech Sat Oct 29 14:43:14 2005| | ||
|LOCATION: sevcech Sat Oct 29 14:43:14 2005| | |-----------------------------------------------------------------------------| | ||
|-----------------------------------------------------------------------------| | | north: -957500 east: -763855.0602047 | | ||
| north: -957500 east: -763855.0602047 | | |REGION south: -1007318 west: -830529 | | ||
|REGION south: -1007318 west: -830529 | | | res: 29.99277544 res: 29.99277544 | | ||
| res: 29.99277544 res: 29.99277544 | | |-----------------------------------------------------------------------------| | ||
|-----------------------------------------------------------------------------| | |MASK:none | | ||
|MASK:none | | |-----------------------------------------------------------------------------| | ||
|-----------------------------------------------------------------------------| | |MAP: (untitled) (ndvi in student) | | ||
|MAP: (untitled) (ndvi in student) | | |-----------------------------------------------------------------------------| | ||
|-----------------------------------------------------------------------------| | | Category Information | cell| % | | | ||
| Category Information | cell| % | | | | #|description | count| cover| hectares| | ||
| #|description | count| cover| hectares| | |-----------------------------------------------------------------------------| | ||
|-----------------------------------------------------------------------------| | | -1--0.925547|from to . . . . . . . . . . . | 37| 0.00| 3.33| | ||
| -1--0.925547|from to . . . . . . . . . . . | 37| 0.00| 3.33| | |-0.553285--0.478832|from to . . . . . . . . . . . | 1| 0.00| 0.09| | ||
|-0.553285--0.478832|from to . . . . . . . . . . . | 1| 0.00| 0.09| | | -0.478832--0.40438|from to . . . . . . . . . . . | 5| 0.00| 0.45| | ||
| -0.478832--0.40438|from to . . . . . . . . . . . | 5| 0.00| 0.45| | | -0.40438--0.329927|from to . . . . . . . . . . . | 114| 0.00| 10.26| | ||
| -0.40438--0.329927|from to . . . . . . . . . . . | 114| 0.00| 10.26| | |-0.329927--0.255474|from to . . . . . . . . . . . | 2258| 0.06| 203.12| | ||
|-0.329927--0.255474|from to . . . . . . . . . . . | 2258| 0.06| 203.12| | |-0.255474--0.181022|from to . . . . . . . . . . . | 13612| 0.38| 1224.49| | ||
|-0.255474--0.181022|from to . . . . . . . . . . . | 13612| 0.38| 1224.49| | |-0.181022--0.106569|from to . . . . . . . . . . . | 11660| 0.32| 1048.89| | ||
|-0.181022--0.106569|from to . . . . . . . . . . . | 11660| 0.32| 1048.89| | |-0.106569--0.032117|from to . . . . . . . . . . . | 70927| 1.96| 6380.36| | ||
|-0.106569--0.032117|from to . . . . . . . . . . . | 70927| 1.96| 6380.36| | | -0.032117-0.042336|from to . . . . . . . . . . . | 88807| 2.46| 7988.78| | ||
| -0.032117-0.042336|from to . . . . . . . . . . . | 88807| 2.46| 7988.78| | | 0.042336-0.116788|from to . . . . . . . . . . . | 171190| 4.74| 15,399.68| | ||
| 0.042336-0.116788|from to . . . . . . . . . . . | 171190| 4.74| 15,399.68| | | 0.116788-0.191241|from to . . . . . . . . . . . | 240714| 6.67| 21,653.83| | ||
| 0.116788-0.191241|from to . . . . . . . . . . . | 240714| 6.67| 21,653.83| | | 0.191241-0.265693|from to . . . . . . . . . . . | 265466| 7.35| 23,880.43| | ||
| 0.191241-0.265693|from to . . . . . . . . . . . | 265466| 7.35| 23,880.43| | | 0.265693-0.340146|from to . . . . . . . . . . . | 342896| 9.50| 30,845.78| | ||
| 0.265693-0.340146|from to . . . . . . . . . . . | 342896| 9.50| 30,845.78| | | 0.340146-0.414599|from to . . . . . . . . . . . | 594617| 16.47| 53,489.76| | ||
| 0.340146-0.414599|from to . . . . . . . . . . . | 594617| 16.47| 53,489.76| | | 0.414599-0.489051|from to . . . . . . . . . . . | 822636| 22.79| 74,001.59| | ||
| 0.414599-0.489051|from to . . . . . . . . . . . | 822636| 22.79| 74,001.59| | | 0.489051-0.563504|from to . . . . . . . . . . . | 752950| 20.86| 67,732.87| | ||
| 0.489051-0.563504|from to . . . . . . . . . . . | 752950| 20.86| 67,732.87| | | 0.563504-0.637956|from to . . . . . . . . . . . | 223519| 6.19| 20,107.02| | ||
| 0.563504-0.637956|from to . . . . . . . . . . . | 223519| 6.19| 20,107.02| | | 0.637956-0.712409|from to . . . . . . . . . . . | 8126| 0.23| 730.99| | ||
| 0.637956-0.712409|from to . . . . . . . . . . . | 8126| 0.23| 730.99| | | 0.712409-0.786861|from to . . . . . . . . . . . | 1| 0.00| 0.09| | ||
| 0.712409-0.786861|from to . . . . . . . . . . . | 1| 0.00| 0.09| | | 0.786861-0.861314|from to . . . . . . . . . . . | 1| 0.00| 0.09| | ||
| 0.786861-0.861314|from to . . . . . . . . . . . | 1| 0.00| 0.09| | |-----------------------------------------------------------------------------| | ||
|-----------------------------------------------------------------------------| | |TOTAL |3609537|100.00|324,701.88| | ||
|TOTAL |3609537|100.00|324,701.88| | +-----------------------------------------------------------------------------+ | ||
+-----------------------------------------------------------------------------+ | |||
Tabulka barev: | Tabulka barev: | ||
-1.00 red | |||
-1.00 red | 0.00 red | ||
0.30 yellow | |||
0.50 green | |||
0.60 0 136 26 | |||
0.87 0 87 16 | |||
#nastavit tabulku barev ze souboru tb-ndvi.txt | |||
#nastavit tabulku barev ze souboru tb-ndvi.txt | # | ||
# | GRASS:~ > r.colors ndvi col=rules < tb-ndvi.txt | ||
GRASS:~ > r.colors ndvi col=rules < tb-ndvi.txt | # | ||
#... a překreslit obsah grafického okna | #... a překreslit obsah grafického okna | ||
# | # | ||
GRASS:~ > d.redraw | GRASS:~ > d.redraw | ||
[[Soubor:ZOD-cv4-ndvi-tb.png|frame|center|Obr č.5: NDVI - vlastní tabulka barev]] | [[Soubor:ZOD-cv4-ndvi-tb.png|frame|center|Obr č.5: NDVI - vlastní tabulka barev]] | ||
Řádek 540: | Řádek 532: | ||
<tt>d.rast.leg</tt>. | <tt>d.rast.leg</tt>. | ||
#současné zobrazení rastrových dat a odpovídající legendy | |||
#současné zobrazení rastrových dat a odpovídající legendy | # | ||
# | GRASS:~ > d.rast.leg map=ndvi | ||
GRASS:~ > d.rast.leg map=ndvi | |||
[[Soubor:ZOD-cv4-ndvi-legenda.png|frame|center|Obr č.6: Zobrazení NDVI a legendy]] | [[Soubor:ZOD-cv4-ndvi-legenda.png|frame|center|Obr č.6: Zobrazení NDVI a legendy]] | ||
Řádek 553: | Řádek 543: | ||
pro ilustraci: | pro ilustraci: | ||
#rozdělit grafické okno na dva rámce | |||
#rozdělit grafické okno na dva rámce | # | ||
GRASS:~ > d.frame f1 at=0,100,0,15;d.frame f2 at=0,100,15,100 | GRASS:~ > d.frame f1 at=0,100,0,15;d.frame f2 at=0,100,15,100 | ||
#v jednom zobrazit ndvi a popisek | # | ||
# | #v jednom zobrazit ndvi a popisek | ||
GRASS:~ > echo "NDVI" | d.text size=5 col=black;d.rast ndvi | # | ||
#a v druhém legendu (pouze pro rozsah [0;1]) | GRASS:~ > echo "NDVI" | d.text size=5 col=black;d.rast ndvi | ||
# | # | ||
GRASS:~ > d.frame -s f1;d.legend ndvi range=0,1 labelnum=7 | #a v druhém legendu (pouze pro rozsah [0;1]) | ||
# | |||
GRASS:~ > d.frame -s f1;d.legend ndvi range=0,1 labelnum=7 | |||
[[Soubor:ZOD-cv4-ndvi-legenda1.png|frame|center|Obr č.7: Zobrazení NDVI a podruhé]] | [[Soubor:ZOD-cv4-ndvi-legenda1.png|frame|center|Obr č.7: Zobrazení NDVI a podruhé]] | ||
Řádek 577: | Řádek 568: | ||
přesnost. | přesnost. | ||
#trik pro r.reclass | |||
#trik pro r.reclass | # | ||
# | GRASS:~ > r.mapcalc 'temp1 = 100 * ndvi' | ||
GRASS:~ > r.mapcalc 'temp1 = 100 * ndvi' | |||
Pro kontrolu vytiskneme rozsah hodnot rastrové vrstvy <u>temp1</u>: | Pro kontrolu vytiskneme rozsah hodnot rastrové vrstvy <u>temp1</u>: | ||
GRASS:~ > r.info temp1 | grep Range | |||
GRASS:~ > r.info temp1 | grep Range | |||
| Range of data: min = -100.000000 max = 86.131387 | |||
| Range of data: min = -100.000000 max = 86.131387 | |||
Připravíme si tzv. ''reklasifikační tabulku'' | Připravíme si tzv. ''reklasifikační tabulku'' | ||
([http:// | ([http://gama.fsv.cvut.cz/~landa/geowikicz_data/zod/cv4/reklas-ndvi.txt reklas-ndvi.txt]), pro bližší informace odkazuji | ||
na manuálové stránky modulu <tt>r.reclass</tt>: | na manuálové stránky modulu <tt>r.reclass</tt>: | ||
-100 thru 5 = 1 bez vegetace, vodni plochy | |||
-100 thru 5 = 1 bez vegetace, vodni plochy | 5 thru 35 = 2 plochy s minimalni vegetaci | ||
35 thru 87 = 3 plochy pokryte vegetaci | |||
Na tomto místě se sluší krátký popis: na každém řádku je uvedeno jedno reklasifikační | Na tomto místě se sluší krátký popis: na každém řádku je uvedeno jedno reklasifikační | ||
Řádek 607: | Řádek 591: | ||
v intervalu od -100 do 5 kategorii (hodnotě) 1, poté následuje textový popisek | v intervalu od -100 do 5 kategorii (hodnotě) 1, poté následuje textový popisek | ||
dané kategorie. Ve svém důsledku bude reklasifikovaná, výstupní rastrová mapa obsahovat | dané kategorie. Ve svém důsledku bude reklasifikovaná, výstupní rastrová mapa obsahovat | ||
celkem tři kategorie (tj. hodnoty). Reklasifikaci provedeme příkazem | celkem tři kategorie (tj. hodnoty). Reklasifikaci provedeme příkazem (předpokládejme, že | ||
máme reklasifikační tabulku předpřipravenu v textovém souboru: | máme reklasifikační tabulku předpřipravenu v textovém souboru): | ||
#reklasifikace rastrových dat | |||
#reklasifikace rastrových dat | # | ||
# | GRASS:~ > r.reclass input=temp1 output=r_ndvi < reklas-ndvi.txt | ||
GRASS:~ > r.reclass input=temp1 output=r_ndvi < reklas-ndvi.txt | |||
Před samotným zobrazením reklasifikovaných dat nastavíme vhodnou tabulku barev | Před samotným zobrazením reklasifikovaných dat nastavíme vhodnou tabulku barev | ||
([http:// | ([http://gama.fsv.cvut.cz/~landa/geowikicz_data/zod/cv4/tb-r_ndvi.txt tb-r_ndvi.txt]): | ||
1 red | |||
1 red | 2 yellow | ||
2 yellow | 3 0 136 26 | ||
3 0 136 26 | |||
#nastavení vhodné tabulky barev | |||
#nastavení vhodné tabulky barev | # | ||
# | GRASS:~ > r.colors r_ndvi col=rules < tb-r_ndvi.txt | ||
GRASS:~ > r.colors r_ndvi col=rules < tb-r_ndvi.txt | # | ||
#zobrazení reklasifikovaných dat | #zobrazení reklasifikovaných dat | ||
GRASS:~ > d.rast r_ndvi | # | ||
GRASS:~ > d.rast.leg r_ndvi | GRASS:~ > d.rast r_ndvi | ||
GRASS:~ > d.rast.leg r_ndvi | |||
[[Soubor:ZOD-cv4-r_ndvi.png|frame|center|Obr č.8: Reklasifikovaný NDVI]] | [[Soubor:ZOD-cv4-r_ndvi.png|frame|center|Obr č.8: Reklasifikovaný NDVI]] | ||
Řádek 646: | Řádek 626: | ||
A ještě krátkou statistiku reklasifikované mapy: | A ještě krátkou statistiku reklasifikované mapy: | ||
GRASS:~ > r.report map=r_ndvi units=c,p,h | |||
GRASS:~ > r.report map=r_ndvi units=c,p,h | |||
+-----------------------------------------------------------------------------+ | +-----------------------------------------------------------------------------+ | ||
| RASTER MAP CATEGORY REPORT | | | RASTER MAP CATEGORY REPORT | | ||
|LOCATION: sevcech Mon Oct 31 01:15:33 2005| | |LOCATION: sevcech Mon Oct 31 01:15:33 2005| | ||
|-----------------------------------------------------------------------------| | |-----------------------------------------------------------------------------| | ||
| north: -957500 east: -763855.0602047 | | | north: -957500 east: -763855.0602047 | | ||
|REGION south: -1007318 west: -830529 | | |REGION south: -1007318 west: -830529 | | ||
| res: 29.99277544 res: 29.99277544 | | | res: 29.99277544 res: 29.99277544 | | ||
|-----------------------------------------------------------------------------| | |-----------------------------------------------------------------------------| | ||
|MASK:none | | |MASK:none | | ||
|-----------------------------------------------------------------------------| | |-----------------------------------------------------------------------------| | ||
|MAP: Reclass of temp1 in student (r_ndvi in student) | | |MAP: Reclass of temp1 in student (r_ndvi in student) | | ||
|-----------------------------------------------------------------------------| | |-----------------------------------------------------------------------------| | ||
| Category Information | cell| % | | | | Category Information | cell| % | | | ||
|#|description | count| cover| hectares| | |#|description | count| cover| hectares| | ||
|-----------------------------------------------------------------------------| | |-----------------------------------------------------------------------------| | ||
|1|bez vegetace, vodni plochy . . . . . . . . . . . | 146375| 3.96| 13,167.41| | |1|bez vegetace, vodni plochy . . . . . . . . . . . | 146375| 3.96| 13,167.41| | ||
|2|plochy s minimalni vegetaci. . . . . . . . . . . | 820825| 22.23| 73,838.67| | |2|plochy s minimalni vegetaci. . . . . . . . . . . | 820825| 22.23| 73,838.67| | ||
|3|plochy pokryte vegetaci. . . . . . . . . . . . . |2642337| 71.56|237,695.81| | |3|plochy pokryte vegetaci. . . . . . . . . . . . . |2642337| 71.56|237,695.81| | ||
|*|no data. . . . . . . . . . . . . . . . . . . . . | 82866| 2.24| 7454.35| | |*|no data. . . . . . . . . . . . . . . . . . . . . | 82866| 2.24| 7454.35| | ||
|-----------------------------------------------------------------------------| | |-----------------------------------------------------------------------------| | ||
|TOTAL |3692403|100.00|332,156.23| | |TOTAL |3692403|100.00|332,156.23| | ||
+-----------------------------------------------------------------------------+ | +-----------------------------------------------------------------------------+ | ||
== praktická úloha: vytvoření mapy s teplotou povrchu == | == praktická úloha: vytvoření mapy s teplotou povrchu == | ||
Řádek 680: | Řádek 659: | ||
v aktuálním mapsetu (student) pod jménem <u>tm6_null</u>: | v aktuálním mapsetu (student) pod jménem <u>tm6_null</u>: | ||
#kopie tm6@PERMANENT -> tm6_null@student | |||
#kopie tm6@PERMANENT -> tm6_null@student | # | ||
# | GRASS:~ > g.copy rast=tm6,tm6_null | ||
GRASS:~ > g.copy rast=tm6,tm6_null | |||
Nyní nahradíme hodnotu "0" hodnotou NULL (žádná data) - nulové hodnoty | Nyní nahradíme hodnotu "0" hodnotou NULL (žádná data) - nulové hodnoty | ||
Řádek 691: | Řádek 668: | ||
příkazem: | příkazem: | ||
#pro jistotu nastavit aktivní region | |||
#pro jistotu nastavit aktivní region | # | ||
# | GRASS:~ > g.region rast=tm6 | ||
GRASS:~ > g.region rast=tm6 | # | ||
#nahrazení hodnoty "0" hodnotou NULL | #nahrazení hodnoty "0" hodnotou NULL | ||
# | # | ||
GRASS:~ > r.null map=tm6_null setnull=0 | GRASS:~ > r.null map=tm6_null setnull=0 | ||
Dle zadaných vzorců (převzato | Dle zadaných vzorců (převzato | ||
z [http://mpa.itc.it/grasstutor/index.phtml Open source GIS: A GRASS GIS Approach]) | z [http://mpa.itc.it/grasstutor/index.phtml Open source GIS: A GRASS GIS Approach]) | ||
provedeme [http:// | provedeme [http://gama.fsv.cvut.cz/~landa/geowikicz_data/zod/cv4/temp-rovnice.txt sérii] výpočtů: | ||
GRASS:~ > r.mapcalc 'tm6rad = 1.238 + (15.6 - 1.238) * tm6_null / 255' | |||
GRASS:~ > r.mapcalc 'tm6rad = 1.238 + (15.6 - 1.238) * tm6_null / 255' | GRASS:~ > r.mapcalc 'temp_kelvin = 1260.56 / (log (607.76 / tm6rad + 1))' | ||
GRASS:~ > r.mapcalc 'temp_kelvin = 1260.56 / (log (607.76 / tm6rad + 1))' | GRASS:~ > r.mapcalc 'temp_celsius = temp_kelvin - 273.15' | ||
GRASS:~ > r.mapcalc 'temp_celsius = temp_kelvin - 273.15' | |||
Před vizualizací výsledné mapy nastavíme vhodnou tabulku barev | Před vizualizací výsledné mapy nastavíme vhodnou tabulku barev | ||
([http:// | ([http://gama.fsv.cvut.cz/~landa/geowikicz_data/zod/cv4/tb-teplota.txt tb-teplota.txt]): | ||
GRASS:~ > r.colors map=temp_celsius col=rules < tb-teplota.txt | |||
GRASS:~ > r.colors map=temp_celsius col=rules < tb-teplota.txt | GRASS:~ > d.rast temp_celsius | ||
GRASS:~ > d.rast temp_celsius | |||
[[Soubor:ZOD-cv4-temp_celsius.png|frame|center|Obr č.10: Teplota povrchu]] | [[Soubor:ZOD-cv4-temp_celsius.png|frame|center|Obr č.10: Teplota povrchu]] | ||
Řádek 726: | Řádek 698: | ||
vznikly během výpočtu <u>temp_celsius</u>: | vznikly během výpočtu <u>temp_celsius</u>: | ||
#odstranění zvolených rastrových mapových vrstev | |||
#odstranění zvolených rastrových mapových vrstev | # | ||
# | GRASS:~ > g.remove rast=tm6_null,tm6rad,temp_kelvin | ||
GRASS:~ > g.remove rast=tm6_null,tm6rad,temp_kelvin | |||
[ [[Zpracování obrazových dat]] ] |
Aktuální verze z 19. 11. 2008, 18:54
Mapová algebra, sčítání, odčítání rastrových dat, podíl obrazu, NDVI
osnova
Úvod do rastrové algebry, výpočty orientované na rastrové buňky - podíl obrazu, normalizovaný diferenční vegetační index (NDVI).
seznam použitých příkazů
- r.mapcalc
- d.rast
- r.info
- r.report
- r.colors
- d.redraw
- d.legend
- d.rast.leg
- r.reclass
- g.copy
- r.null
- g.remove
rastrová algebra
Rastrová (mapová) algebra je poměrně komplexní téma. Vzhledem k rozsahu problematiky si v rámci tohoto kurzu nastíníme pouze její základy. Naše znalosti v této oblasti zpracování obrazu budeme postupně během následujících cvičení dále prohlubovat. Rastrovou algebru má v GRASSu ve své moci komplexní modul r.mapcalc. Jde o "mocný" modul s rozsáhlou manuálovou stránkou (viz g.manual r.mapcalc). K tématu rastrové algebry je dostupný i speciální tutoriál. Rastrové algebře je věnována kapitola 6.16 v GRASS příručce a kapitola 15 v GIS GRASS 6.0 - Praktická rukověť začínajících uživatelů. Pro lepší pochopení probírané látky, doporučuji obě dvě kapitoly alespoň zběžně prostudovat.
Chování modulu r.mapcalc je celkem standardní. Lze ho ovládat ve dvou různým módech - tzv. interaktivním, kdy spustíte modul bez parametrů a zadáváte jednotlivé výrazy, pro ukončení slouží klíčové slovo "end". Tedy podobně jako při zadávání pravidel pro vlastní tabulku barev v modulu r.colors.
#interaktivní mód modulu r.mapcalc # GRASS:~ > r.mapcalc Enter expressions, "end" when done. mapcalc> p12 = tm1 / tm2 mapcalc> p23 = tm2 / tm3 mapcalc> end
Výhodnější je však předat dané výrazy modulu jako parametr (uzavřít do uvozovek!), nebo je uložit do textového souboru (například při opětovném použití) a přesměrovat do modulu přes standardní vstup (podobně jako tabulku barev do r.colors). Např. (podily.txt):
#zadání rovnice jako parametr # GRASS:~ > r.mapcalc 'p12 = tm1 / tm2' GRASS:~ > r.mapcalc 'p23 = tm2 / tm3' # #rovnice uložené v textovém souboru # GRASS:~ > r.mapcalc < podily.txt
My bude upřednostňovat způsob zadávání výrazů jako parametru.
Poznámka: Mějte na paměti, že i tento modul provádí výpočet pouze v aktuální regionu!
Obecně výraz pro r.mapcalc vypadá následovně:
nova_mapa = vyraz (mapa1, mapa2, ...)
kde vyraz může zahrnovat názvy rastrových mapových vrstev, celočíselné či reálné konstanty a různé funkce.
operátor | popis | typ | priorita |
---|---|---|---|
^ | exponent | aritmetický | 5 |
% | modulo (zbytek po celočíselném dělení) | aritmetický | 4 |
/ | dělení | aritmetický | 4 |
* | násobení | aritmetický | 4 |
+ | sčítání | aritmetický | 3 |
- | odečítání | aritmetický | 3 |
== | rovnost | logický | 2 |
!= | nerovnost | logický | 2 |
> | větší než | logický | 2 |
< | menší než | logický | 2 |
>= | větší rovno než | logický | 2 |
<= | menší rovno než | logický | 2 |
&& | AND (a) | logický | 1 |
|| | OR (nebo) | logický | 1 |
# | rozděluje tabulku barev | aritmetický | - |
funkce | význam | mat. zápis | typ |
---|---|---|---|
abs(x) | absolutní hodnota z x | |x| | |
atan(x) | arkustangens z x (výsledek ve stupních) | arctan(x) | F |
atan(x,y) | arkustangens z x/y (výsledek ve stupních) | arctan(x/y) | F |
cos(x) | kosinus z x (x ve stupních) | cos(x) | F |
double(x) | převede x na číslo s plovoucí desetinnou čárkou (výstupní mapa DCELL) | F | |
eval([x,y,...,]z) | vyhodnotí výrazy x a y a výsledek uloží do z | ||
exp(x) | exponenciální funkce z x | ex | F |
exp(x,y) | mocnina y z x | xy | F |
float(x) | převede x na číslo s plovoucí desetinnou čárkou (výstupní mapa FCELL) | F | |
if | podmíněný příkaz | ||
if(x) | 1, když x není 0, jinak 0 | ||
if(x,a) | a, když x není 0, jinak 0 | ||
if(x,a,b) | a, když x není 0, jinak b | ||
if(x,a,b,c) | a, když x > 0, b když x rovno 0, c když x < 0 | ||
int(x) | převede x na celé číslo, zbytek odřízne | ||
isnull() | 1, když x rovno NULL (no data) | F | |
log(x) | přirozený logaritmus z x | ln(x) | F |
log(x,b) | logaritmus z x při základu b | logb(x) | F |
max(x,y [,z...]) | vrátí největší ze zadaných hodnot | * | |
median(x,y[,z...]) | vrátí prostřední hodnotu (medián) ze zadaných hodnot | * | |
min(x,y[,z...]) | vrátí nejmenší ze zadaných hodnot | * | |
mode(x,y[,z...]) | vrátí nejčetnější hodnotu (modus) ze zadaných hodnot | * | |
not(x) | 1 pokud x je nula, jinak 0 | ||
rand(low,high) | vrátí náhodnou číselnou hodnotu z intervalu [low; high] | ||
round(x) | zaokrouhlí x na nejbližší celočíselnou hodnotu | I | |
sin(x) | sinus z x (x ve stupních) | sin(x) | F |
sqrt(x) | odmocnina z x | sqrt(x) | F |
tan(x) | tangens z x (x ve stupních) | tan(x) | F |
Vysvětlivky:
* | ... | Vrací reálné číslo, pokud je x reálné číslo. |
F | ... | Výsledek je vždy reálné číslo. |
I | ... | Výsledek je celé číslo. |
proměnná | význam |
---|---|
x() | proměnná pro aktuální souřadnici x |
y() | proměnná pro aktuální souřadnici y |
col() | proměnná pro aktuální číslo sloupce |
row() | proměnná pro aktuální číslo řádku |
nsres() | proměnná pro aktuální rozlišení sever-jih |
ewres() | proměnná pro aktuální rozlišení východ-západ |
V této části budeme hovořit pouze o výpočtech pixelově orientovaných (hodnota výsledné buňky vzniká na základě hodnot jednotlivých pixelů vstupní mapy/map). V některém z příštích kurzů narazíme na výpočty, které se provádí v rámci množiny pixelů pokryté tzv. kernelem - plovoucím oknem.
Vybereme-li dostatečný detail, můžeme velmi snadno výsledky zvolených operací "překontrolovat".
Několik drobných poznámek k chování modulu:
- Typ výsledné mapy (celočíselná; hodnoty s plovoucí desetinnou čárkou) závisí na typu vstupních map a datovém typu konstant uvedených ve výrazu. Výsledkem dělení dvou celočíselných map je tedy opět mapa celočíselná!
- Vstupují-li do výpočtu rastrové buňky s hodnotou NULL (žádná data) je výsledkem opět hodnota NULL.
Asi nemá význam se dále potácet na poli pusté teorie, přejděme tedy raději k praktickým úlohám. Modul r.mapcalc budeme v dalších kurzech hojně používat, naši pozornost si rozhodně zaslouží...
podíl obrazu
Podíl obrazu je základní a velmi jednoduchou algebraickou metodou používanou pro detekci příznaků, zvýraznění obrazu, výpočet vegetačních indexů a řadu dalších úloh. Prakticky se podíl obrazu provádí pomocí komplexního modulu pro mapovou algebru - r.mapcalc. Velmi důležité je doplnit rovnici o násobení konstantou "1.0". Dělíme totiž celočíselná rastrová data. Výsledek je tedy, jak již víme, opět celočíselný! Trik s vynásobením hodnotou "1.0" má za následek vytvoření výstupní rastrové vrstvy obsahující hodnoty s plovoucí desetinnou čárkou (DCELL).
Jako příklad můžeme vypočítat podíl sedmého a čtvrtého kanálu družicové scény LandSat-TM5:
#jistota... # GRASS:~ > g.region rast=tm4 # #podíl obrazových kanálů sedm a čtyři (výsledná mapa DCELL) # GRASS:~ > r.mapcalc 'p74 = 1.0 * tm7 / tm4' GRASS:~ > d.rast p74
Praktickou úlohou na podíl obrazu je výpočet normalizovaného diferenčního vegetačního indexu (NDVI). Příkaz pro tvorbu mapy s NDVI vypadá následovně:
#výpočet NDVI # GRASS:~ > r.mapcalc 'ndvi = 1.0 * (tm4 - tm3) / (tm4 + tm3)'
Dominuje-li DH v blízkém infračerveném oboru (tm4) nad hodnotou ve viditelném (červeném) spektru (tm3), hodnota NDVI je kladná a ukazuje na oblast s bohatou vegetací. Hodnoty NDVI kolem nuly reprezentují půdu bez výraznější vegetace, nula potom vodní plochy - pouze teoreticky.
Pro lepší vizuální interpretaci obrazu zvolíme vlastní tabulku barev (tb-ndvi.txt), nejprve zjistíme rozsah dat (r.info), pro názornější představu můžeme vytisknout statistiku dat (r.report):
#rozsah hodnot # GRASS:~ > r.info ndvi +----------------------------------------------------------------------------+ | Layer: ndvi Date: Sat Oct 29 14:42:15 2005 | | Mapset: student Login of Creator: martin | | Location: sevcech | | DataBase: /home/martin/grassdata | | Title: ( ndvi ) | |----------------------------------------------------------------------------| | | | Type of Map: cell Number of Categories: 255 | | Data Type: DCELL | | Rows: 1661 | | Columns: 2223 | | Total Cells: 3692403 | | Projection: krovak (zone 0) | | N: -957500 S: -1007318 Res: 29.99277544 | | E: -763855.0602047 W: -830529 Res: 29.99277544 | | Range of data: min = -1.000000 max = 0.861314 | | | | Data Source: | | | | | | | | Data Description: | | generated by r.mapcalc | | | | Comments: | | 1.000000 * (tm4 - tm3) / (tm4 + tm3) | | | +----------------------------------------------------------------------------+
Rastrové hodnoty se tedy pohybují v intervalu od [-1.0; +0.87]. Všimněme si navíc typu dat - DCELL (rastrová data s plovoucí desetinnou čárkou). V komentáři (Comments) je navíc uveden výraz, na základě něhož tato data vznikla.
#vytisknout statistiku (počet buněk, procenta, plocha v hektarech); přerozdělit na 25 intervalů # GRASS:~ > r.report -n map=ndvi units=c,p,h nstep=25 +-----------------------------------------------------------------------------+ | RASTER MAP CATEGORY REPORT | |LOCATION: sevcech Sat Oct 29 14:43:14 2005| |-----------------------------------------------------------------------------| | north: -957500 east: -763855.0602047 | |REGION south: -1007318 west: -830529 | | res: 29.99277544 res: 29.99277544 | |-----------------------------------------------------------------------------| |MASK:none | |-----------------------------------------------------------------------------| |MAP: (untitled) (ndvi in student) | |-----------------------------------------------------------------------------| | Category Information | cell| % | | | #|description | count| cover| hectares| |-----------------------------------------------------------------------------| | -1--0.925547|from to . . . . . . . . . . . | 37| 0.00| 3.33| |-0.553285--0.478832|from to . . . . . . . . . . . | 1| 0.00| 0.09| | -0.478832--0.40438|from to . . . . . . . . . . . | 5| 0.00| 0.45| | -0.40438--0.329927|from to . . . . . . . . . . . | 114| 0.00| 10.26| |-0.329927--0.255474|from to . . . . . . . . . . . | 2258| 0.06| 203.12| |-0.255474--0.181022|from to . . . . . . . . . . . | 13612| 0.38| 1224.49| |-0.181022--0.106569|from to . . . . . . . . . . . | 11660| 0.32| 1048.89| |-0.106569--0.032117|from to . . . . . . . . . . . | 70927| 1.96| 6380.36| | -0.032117-0.042336|from to . . . . . . . . . . . | 88807| 2.46| 7988.78| | 0.042336-0.116788|from to . . . . . . . . . . . | 171190| 4.74| 15,399.68| | 0.116788-0.191241|from to . . . . . . . . . . . | 240714| 6.67| 21,653.83| | 0.191241-0.265693|from to . . . . . . . . . . . | 265466| 7.35| 23,880.43| | 0.265693-0.340146|from to . . . . . . . . . . . | 342896| 9.50| 30,845.78| | 0.340146-0.414599|from to . . . . . . . . . . . | 594617| 16.47| 53,489.76| | 0.414599-0.489051|from to . . . . . . . . . . . | 822636| 22.79| 74,001.59| | 0.489051-0.563504|from to . . . . . . . . . . . | 752950| 20.86| 67,732.87| | 0.563504-0.637956|from to . . . . . . . . . . . | 223519| 6.19| 20,107.02| | 0.637956-0.712409|from to . . . . . . . . . . . | 8126| 0.23| 730.99| | 0.712409-0.786861|from to . . . . . . . . . . . | 1| 0.00| 0.09| | 0.786861-0.861314|from to . . . . . . . . . . . | 1| 0.00| 0.09| |-----------------------------------------------------------------------------| |TOTAL |3609537|100.00|324,701.88| +-----------------------------------------------------------------------------+
Tabulka barev:
-1.00 red 0.00 red 0.30 yellow 0.50 green 0.60 0 136 26 0.87 0 87 16
#nastavit tabulku barev ze souboru tb-ndvi.txt # GRASS:~ > r.colors ndvi col=rules < tb-ndvi.txt # #... a překreslit obsah grafického okna # GRASS:~ > d.redraw
Na tomto místě si ukážeme, jak zobrazit legendu k rastrové vrstvě. Pro zobrazení samostatné legendy slouží modul d.legend, současné zobrazení legendy a rastrové mapové vrstvy nabízí (za cenu omezené možnosti ovlivnit podobu legendy) modul d.rast.leg.
#současné zobrazení rastrových dat a odpovídající legendy # GRASS:~ > d.rast.leg map=ndvi
Poznámka: Jelikož d.rast.leg rozdělí GRASS monitor na dvě různé části (rámce), je nutné jej posléze vyčistit příkazem d.erase -f.
Požadujeme-li více flexibilní výstup, je to pěkná hromádka příkazů, pouze pro ilustraci:
#rozdělit grafické okno na dva rámce # GRASS:~ > d.frame f1 at=0,100,0,15;d.frame f2 at=0,100,15,100 # #v jednom zobrazit ndvi a popisek # GRASS:~ > echo "NDVI" | d.text size=5 col=black;d.rast ndvi # #a v druhém legendu (pouze pro rozsah [0;1]) # GRASS:~ > d.frame -s f1;d.legend ndvi range=0,1 labelnum=7
reklasifikace rastrových dat
Chceme-li na základě NDVI provést primitivní klasifikaci na tři základní třídy (bez vegetace, vodní plochy / plochy s minimální vegetací / plochy pokryté vegetací), postupujeme takto:
Protože modul určený pro reklasifikaci rastrových dat r.reclass nepodporuje data s plovoucí desetinnou čárkou, musíme si pomoci malým trikem - vytvoříme pomocnou vrstvu, která vznikne přenásobením NDVI konstantou "100". Zachováme si tak dostačující přesnost.
#trik pro r.reclass # GRASS:~ > r.mapcalc 'temp1 = 100 * ndvi'
Pro kontrolu vytiskneme rozsah hodnot rastrové vrstvy temp1:
GRASS:~ > r.info temp1 | grep Range | Range of data: min = -100.000000 max = 86.131387
Připravíme si tzv. reklasifikační tabulku (reklas-ndvi.txt), pro bližší informace odkazuji na manuálové stránky modulu r.reclass:
-100 thru 5 = 1 bez vegetace, vodni plochy 5 thru 35 = 2 plochy s minimalni vegetaci 35 thru 87 = 3 plochy pokryte vegetaci
Na tomto místě se sluší krátký popis: na každém řádku je uvedeno jedno reklasifikační pravidlo, první řádek například "říká": proveď přiřazení hodnot vstupních dat v intervalu od -100 do 5 kategorii (hodnotě) 1, poté následuje textový popisek dané kategorie. Ve svém důsledku bude reklasifikovaná, výstupní rastrová mapa obsahovat celkem tři kategorie (tj. hodnoty). Reklasifikaci provedeme příkazem (předpokládejme, že máme reklasifikační tabulku předpřipravenu v textovém souboru):
#reklasifikace rastrových dat # GRASS:~ > r.reclass input=temp1 output=r_ndvi < reklas-ndvi.txt
Před samotným zobrazením reklasifikovaných dat nastavíme vhodnou tabulku barev (tb-r_ndvi.txt):
1 red 2 yellow 3 0 136 26
#nastavení vhodné tabulky barev # GRASS:~ > r.colors r_ndvi col=rules < tb-r_ndvi.txt # #zobrazení reklasifikovaných dat # GRASS:~ > d.rast r_ndvi GRASS:~ > d.rast.leg r_ndvi
Důležité je poznamenat, že reklasifikační modul nevytváří ve skutečnosti plnohodnotnou rastrovou vrstvu (šlo by o zbytečné plýtvání diskového prostoru), pouze zapíše zvolená reklasifikační pravidla. Tento fakt je dobré mít na paměti, například to vede k tomu, že nelze smazat podkladovou vrstvu, která sloužila jako vstup do reklasifikace před samotnou reklasifikovanou mapou.
A ještě krátkou statistiku reklasifikované mapy:
GRASS:~ > r.report map=r_ndvi units=c,p,h +-----------------------------------------------------------------------------+ | RASTER MAP CATEGORY REPORT | |LOCATION: sevcech Mon Oct 31 01:15:33 2005| |-----------------------------------------------------------------------------| | north: -957500 east: -763855.0602047 | |REGION south: -1007318 west: -830529 | | res: 29.99277544 res: 29.99277544 | |-----------------------------------------------------------------------------| |MASK:none | |-----------------------------------------------------------------------------| |MAP: Reclass of temp1 in student (r_ndvi in student) | |-----------------------------------------------------------------------------| | Category Information | cell| % | | |#|description | count| cover| hectares| |-----------------------------------------------------------------------------| |1|bez vegetace, vodni plochy . . . . . . . . . . . | 146375| 3.96| 13,167.41| |2|plochy s minimalni vegetaci. . . . . . . . . . . | 820825| 22.23| 73,838.67| |3|plochy pokryte vegetaci. . . . . . . . . . . . . |2642337| 71.56|237,695.81| |*|no data. . . . . . . . . . . . . . . . . . . . . | 82866| 2.24| 7454.35| |-----------------------------------------------------------------------------| |TOTAL |3692403|100.00|332,156.23| +-----------------------------------------------------------------------------+
praktická úloha: vytvoření mapy s teplotou povrchu
Před samotným výpočtem vytvoříme kopii šestého kanálu LandSat-TM5 (seznámíme se tedy s novým modulem - g.copy určeným pro kopírování jednotlivých mapových vrstev) - tj. kopii tm6 z mapsetu PERMANENT, kterou uložíme v aktuálním mapsetu (student) pod jménem tm6_null:
#kopie tm6@PERMANENT -> tm6_null@student # GRASS:~ > g.copy rast=tm6,tm6_null
Nyní nahradíme hodnotu "0" hodnotou NULL (žádná data) - nulové hodnoty obsahuje levý dolní roh snímku, který nenese žádná data. Nahradíme-li tuto hodnotu hodnotou NULL, nebude tato část snímku vstupovat do výpočtu. Nahrazení hodnotou NULL provedeme příkazem:
#pro jistotu nastavit aktivní region # GRASS:~ > g.region rast=tm6 # #nahrazení hodnoty "0" hodnotou NULL # GRASS:~ > r.null map=tm6_null setnull=0
Dle zadaných vzorců (převzato z Open source GIS: A GRASS GIS Approach) provedeme sérii výpočtů:
GRASS:~ > r.mapcalc 'tm6rad = 1.238 + (15.6 - 1.238) * tm6_null / 255' GRASS:~ > r.mapcalc 'temp_kelvin = 1260.56 / (log (607.76 / tm6rad + 1))' GRASS:~ > r.mapcalc 'temp_celsius = temp_kelvin - 273.15'
Před vizualizací výsledné mapy nastavíme vhodnou tabulku barev (tb-teplota.txt):
GRASS:~ > r.colors map=temp_celsius col=rules < tb-teplota.txt GRASS:~ > d.rast temp_celsius
Zde se hodí představit další všeobecný modul a to g.remove. Jak jeho název napovídá, slouží pro mazaní jednotlivých mapových vrstev. Nás v podstatě budou zajímat pouze parametry rast pro rastrová a vect pro vektorová data. Nic nám tedy nebrání odstranit již nepotřebné rastrové vrstvy, které vznikly během výpočtu temp_celsius:
#odstranění zvolených rastrových mapových vrstev # GRASS:~ > g.remove rast=tm6_null,tm6rad,temp_kelvin