153YZOD Zpracování obrazových dat 2006 - 4. cvičení: Porovnání verzí

Z GeoWikiCZ
Přejít na: navigace, hledání
(rastrová algebra: - slouceni obrazku (obr c.1))
m (upozorneni)
 
(Není zobrazeno 20 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/>&nbsp;</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__
 
== osnova ==
 
== osnova ==
  
Řádek 8: Řádek 15:
 
== seznam použitých příkazů ==
 
== seznam použitých příkazů ==
  
* <tt>[http://grass.itc.it/grass60/manuals/html60_user/r.mapcalc.html r.mapcalc]</tt>
+
* [http://grass.itc.it/grass60/manuals/html60_user/r.mapcalc.html r.mapcalc]
* <tt>[http://grass.itc.it/grass60/manuals/html60_user/d.rast.html d.rast]</tt>
+
* [http://grass.itc.it/grass60/manuals/html60_user/d.rast.html d.rast]
* <tt>[http://grass.itc.it/grass60/manuals/html60_user/r.info.html r.info]</tt>
+
* [http://grass.itc.it/grass60/manuals/html60_user/r.info.html r.info]
* <tt>[http://grass.itc.it/grass60/manuals/html60_user/r.report.html r.report]</tt>
+
* [http://grass.itc.it/grass60/manuals/html60_user/r.report.html r.report]
* <tt>[http://grass.itc.it/grass60/manuals/html60_user/r.colors.html r.colors]</tt>
+
* [http://grass.itc.it/grass60/manuals/html60_user/r.colors.html r.colors]
* <tt>[http://grass.itc.it/grass60/manuals/html60_user/d.redraw.html d.redraw]</tt>
+
* [http://grass.itc.it/grass60/manuals/html60_user/d.redraw.html d.redraw]
* <tt>[http://grass.itc.it/grass60/manuals/html60_user/d.legend.html d.legend]</tt>
+
* [http://grass.itc.it/grass60/manuals/html60_user/d.legend.html d.legend]
* <tt>[http://grass.itc.it/grass60/manuals/html60_user/d.rast.leg.html d.rast.leg]</tt>
+
* [http://grass.itc.it/grass60/manuals/html60_user/d.rast.leg.html d.rast.leg]
* <tt>[http://grass.itc.it/grass60/manuals/html60_user/r.reclass.html r.reclass]</tt>
+
* [http://grass.itc.it/grass60/manuals/html60_user/r.reclass.html r.reclass]
* <tt>[http://grass.itc.it/grass60/manuals/html60_user/g.copy.html g.copy]</tt>
+
* [http://grass.itc.it/grass60/manuals/html60_user/g.copy.html g.copy]
* <tt>[http://grass.itc.it/grass60/manuals/html60_user/r.null.html r.null]</tt>
+
* [http://grass.itc.it/grass60/manuals/html60_user/r.null.html r.null]
* <tt>[http://grass.itc.it/grass60/manuals/html60_user/g.remove.html g.remove]</tt>
+
* [http://grass.itc.it/grass60/manuals/html60_user/g.remove.html g.remove]
  
 
== rastrová algebra ==
 
== rastrová algebra ==
Řádek 42: Řádek 49:
 
pravidel pro vlastní tabulku barev v&nbsp;modulu <tt>r.colors</tt>.
 
pravidel pro vlastní tabulku barev v&nbsp;modulu <tt>r.colors</tt>.
  
<pre>
+
#interaktivní mód modulu r.mapcalc
#interaktivní mód modulu r.mapcalc
+
#
#
+
GRASS:~&nbsp;>&nbsp;r.mapcalc
GRASS:~&nbsp;>&nbsp;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
</pre>
 
  
 
Výhodnější je však předat dané výrazy modulu jako parametr (uzavřít do&nbsp;<u>uvozovek</u>!), nebo
 
Výhodnější je však předat dané výrazy modulu jako parametr (uzavřít do&nbsp;<u>uvozovek</u>!), nebo
 
je uložit do&nbsp;textového souboru (například při opětovném použití)
 
je uložit do&nbsp;textového souboru (například při opětovném použití)
 
a&nbsp;přesměrovat do&nbsp;modulu přes standardní vstup (podobně jako tabulku barev
 
a&nbsp;přesměrovat do&nbsp;modulu přes standardní vstup (podobně jako tabulku barev
do&nbsp;<tt>r.colors</tt>). Např.&nbsp;([http://grass.fsv.cvut.cz/~martin/vyuka/zzoz/cv4/podily.txt podily.txt]):
+
do&nbsp;<tt>r.colors</tt>). Např.&nbsp;([http://gama.fsv.cvut.cz/~landa/geowikicz_data/zod/cv4/podily.txt podily.txt]):
  
<pre>
+
#zadání rovnice jako parametr
#zadání rovnice jako parametr
+
#
#
+
GRASS:~&nbsp;>&nbsp;r.mapcalc 'p12 = tm1 / tm2'
GRASS:~&nbsp;>&nbsp;r.mapcalc 'p12 = tm1 / tm2'
+
GRASS:~&nbsp;>&nbsp;r.mapcalc 'p23 = tm2 / tm3'
GRASS:~&nbsp;>&nbsp;r.mapcalc 'p23 = tm2 / tm3'
+
#
#rovnice uložené v textovém souboru
+
#rovnice uložené v textovém souboru
#
+
#
GRASS:~&nbsp;>&nbsp;r.mapcalc < podily.txt
+
GRASS:~&nbsp;>&nbsp;r.mapcalc < podily.txt
</pre>
 
  
 
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 74: Řádek 79:
 
Obecně výraz pro <tt>r.mapcalc</tt> vypadá následovně:
 
Obecně výraz pro <tt>r.mapcalc</tt> vypadá následovně:
  
<pre>
 
 
  nova_mapa = vyraz (mapa1, mapa2, ...)
 
  nova_mapa = vyraz (mapa1, mapa2, ...)
</pre>
 
  
 
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 366: Řádek 369:
 
detail, můžeme velmi snadno výsledky zvolených operací "překontrolovat".
 
detail, můžeme velmi snadno výsledky zvolených operací "překontrolovat".
  
[[Soubor:ZOD-cv4-test-tm1_tm2.png|frame|center|Obr č.1 Vstupní kanál tm1 a tm2]]
+
[[Soubor:ZOD-cv4-test-tm1_tm2.png|frame|center|Obr č.1 Vstupní kanály tm1 a tm2]]
  
[[Soubor:ZOD-cv4/test-operace12.png|frame|center|Obr č.2: (a) tm1 + tm2; (b) tm1 - tm2; (c) tm1 * tm2; (d) tm1 / tm2; (e) 1.0 * tm1 / tm2]]
+
[[Soubor:ZOD-cv4-test-operace12.png|frame|center|Obr č.2: (a) tm1 + tm2; (b) tm1 - tm2; (c) tm1 * tm2; (d) tm1 / tm2; (e) 1.0 * tm1 / tm2]]
  
 
Několik drobných poznámek k&nbsp;chování modulu:
 
Několik drobných poznámek k&nbsp;chování modulu:
Řádek 379: Řádek 382:
 
naši pozornost si rozhodně zaslouží...
 
naši pozornost si rozhodně zaslouží...
  
== podíl obrazů ==
+
== podíl obrazu ==
  
 
Podíl obrazu je základní a&nbsp;velmi jednoduchou algebraickou metodou používanou pro detekci příznaků,
 
Podíl obrazu je základní a&nbsp;velmi jednoduchou algebraickou metodou používanou pro detekci příznaků,
Řádek 390: Řádek 393:
 
Jako příklad můžeme vypočítat podíl sedmého a&nbsp;čtvrtého kanálu družicové scény LandSat-TM5:
 
Jako příklad můžeme vypočítat podíl sedmého a&nbsp;čtvrtého kanálu družicové scény LandSat-TM5:
  
<pre>
+
#jistota...
#jistota...
+
#
#
+
GRASS:~&nbsp;>&nbsp;g.region rast=tm4
GRASS:~&nbsp;>&nbsp;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:~&nbsp;>&nbsp;r.mapcalc 'p74 = 1.0 * tm7 / tm4'
+
GRASS:~&nbsp;>&nbsp;r.mapcalc 'p74 = 1.0 * tm7 / tm4'
GRASS:~&nbsp;>&nbsp;d.rast p74
+
GRASS:~&nbsp;>&nbsp;d.rast p74
</pre>
 
  
 
[[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 405: Řádek 407:
 
Příkaz pro tvorbu mapy s&nbsp;NDVI vypadá následovně:
 
Příkaz pro tvorbu mapy s&nbsp;NDVI vypadá následovně:
  
<pre>
+
#výpočet NDVI
#výpočet NDVI
+
#
#
+
GRASS:~&nbsp;>&nbsp;r.mapcalc 'ndvi = 1.0 * (tm4 - tm3) / (tm4 + tm3)'
GRASS:~&nbsp;>&nbsp;r.mapcalc 'ndvi = 1.0 * (tm4 - tm3) / (tm4 + tm3)'
 
</pre>
 
  
 
Dominuje-li DH v&nbsp;blízkém infračerveném oboru (<u>tm4</u>) nad hodnotou
 
Dominuje-li DH v&nbsp;blízkém infračerveném oboru (<u>tm4</u>) nad hodnotou
Řádek 421: Řá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://grass.fsv.cvut.cz/~martin/vyuka/zzoz/cv4/tb-ndvi.txt tb-ndvi.txt]), nejprve zjistíme rozsah dat
+
([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>):
  
<pre>
+
#rozsah hodnot
#rozsah hodnot
+
#
#
+
GRASS:~&nbsp;>&nbsp;r.info ndvi
GRASS:~&nbsp;>&nbsp;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)                                    |
+
|                                                                            |
|                                                                            |
+
+----------------------------------------------------------------------------+
+----------------------------------------------------------------------------+
 
</pre>
 
  
 
Rastrové hodnoty se tedy pohybují v&nbsp;intervalu od [-1.0;&nbsp;+0.87]. Všimněme si navíc
 
Rastrové hodnoty se tedy pohybují v&nbsp;intervalu od [-1.0;&nbsp;+0.87]. Všimněme si navíc
typu dat - DCELL (rastrová data s&nbsp;plovoucí desetinnou čárkou). V&nbsp;komentáři (Comments)
+
typu dat - DCELL (rastrová data s&nbsp;plovoucí desetinnou čárkou). V&nbsp;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.
  
<pre>
+
#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:~&nbsp;>&nbsp;r.report -n map=ndvi units=c,p,h nstep=25
GRASS:~&nbsp;>&nbsp;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|
+
+-----------------------------------------------------------------------------+
+-----------------------------------------------------------------------------+
 
</pre>
 
  
 
Tabulka barev:
 
Tabulka barev:
  
<pre>
+
-1.00  red
-1.00  red
+
  0.00  red
0.00  red
+
  0.30  yellow
0.30  yellow
+
  0.50  green
0.50  green
+
  0.60  0 136 26
0.60  0 136 26
+
  0.87  0  87 16
0.87  0  87 16
 
</pre>
 
  
<pre>
+
#nastavit tabulku barev ze souboru tb-ndvi.txt
#nastavit tabulku barev ze souboru tb-ndvi.txt
+
#
#
+
GRASS:~&nbsp;>&nbsp;r.colors ndvi col=rules < tb-ndvi.txt
GRASS:~&nbsp;>&nbsp;r.colors ndvi col=rules < tb-ndvi.txt
+
#
#... a překreslit obsah grafického okna
+
#... a překreslit obsah grafického okna
#
+
#
GRASS:~&nbsp;>&nbsp;d.redraw
+
GRASS:~&nbsp;>&nbsp;d.redraw
</pre>
 
  
 
[[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 539: Řádek 532:
 
<tt>d.rast.leg</tt>.
 
<tt>d.rast.leg</tt>.
  
<pre>
+
#současné zobrazení rastrových dat a odpovídající legendy
#současné zobrazení rastrových dat a odpovídající legendy
+
#
#
+
GRASS:~&nbsp;>&nbsp;d.rast.leg map=ndvi
GRASS:~&nbsp;>&nbsp;d.rast.leg map=ndvi
 
</pre>
 
  
 
[[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 552: Řádek 543:
 
pro ilustraci:
 
pro ilustraci:
  
<pre>
+
#rozdělit grafické okno na dva rámce
#rozdělit grafické okno na dva rámce
+
#
GRASS:~&nbsp;>&nbsp;d.frame f1 at=0,100,0,15;d.frame f2 at=0,100,15,100</tr>
+
GRASS:~&nbsp;>&nbsp;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:~&nbsp;>&nbsp;echo "NDVI" | d.text size=5 col=black;d.rast ndvi
+
#
#a v druhém legendu (pouze pro rozsah [0;1])
+
GRASS:~&nbsp;>&nbsp;echo "NDVI" | d.text size=5 col=black;d.rast ndvi
#
+
#
GRASS:~&nbsp;>&nbsp;d.frame -s f1;d.legend ndvi range=0,1 labelnum=7
+
#a v druhém legendu (pouze pro rozsah [0;1])
</pre>
+
#
 +
GRASS:~&nbsp;>&nbsp;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 576: Řádek 568:
 
přesnost.
 
přesnost.
  
<pre>
+
#trik pro r.reclass
#trik pro r.reclass
+
#
#
+
GRASS:~&nbsp;>&nbsp;r.mapcalc 'temp1 = 100 * ndvi'
GRASS:~&nbsp;>&nbsp;r.mapcalc 'temp1 = 100 * ndvi'
 
</pre>
 
  
 
Pro kontrolu vytiskneme rozsah hodnot rastrové vrstvy <u>temp1</u>:
 
Pro kontrolu vytiskneme rozsah hodnot rastrové vrstvy <u>temp1</u>:
  
<pre>
+
GRASS:~&nbsp;>&nbsp;r.info temp1 | grep Range
GRASS:~&nbsp;>&nbsp;r.info temp1 | grep Range
+
</pre>
+
 
+
|  Range of data:    min =  -100.000000 max = 86.131387
<pre>
 
|  Range of data:    min =  -100.000000 max = 86.131387
 
</pre>
 
  
 
Připravíme si tzv.&nbsp;''reklasifikační tabulku''
 
Připravíme si tzv.&nbsp;''reklasifikační tabulku''
([http://grass.fsv.cvut.cz/~martin/vyuka/zzoz/cv4/reklas-ndvi.txt reklas-ndvi.txt]), pro bližší informace odkazuji
+
([http://gama.fsv.cvut.cz/~landa/geowikicz_data/zod/cv4/reklas-ndvi.txt reklas-ndvi.txt]), pro bližší informace odkazuji
 
na&nbsp;manuálové stránky modulu <tt>r.reclass</tt>:
 
na&nbsp;manuálové stránky modulu <tt>r.reclass</tt>:
  
<pre>
+
-100 thru 5  = 1 bez vegetace, vodni plochy
-100 thru 5  = 1 bez vegetace, vodni plochy
+
  5  thru 35  = 2 plochy s minimalni vegetaci
5  thru 35  = 2 plochy s minimalni vegetaci
+
  35  thru 87  = 3 plochy pokryte vegetaci
35  thru 87  = 3 plochy pokryte vegetaci
 
</pre>
 
  
 
Na tomto místě se sluší krátký popis: na&nbsp;každém řádku je uvedeno jedno reklasifikační
 
Na tomto místě se sluší krátký popis: na&nbsp;každém řádku je uvedeno jedno reklasifikační
Řádek 606: Řádek 591:
 
v&nbsp;intervalu od -100 do 5 kategorii (hodnotě)&nbsp;1, poté následuje textový popisek
 
v&nbsp;intervalu od -100 do 5 kategorii (hodnotě)&nbsp;1, poté následuje textový popisek
 
dané kategorie. Ve&nbsp;svém důsledku bude reklasifikovaná, výstupní rastrová mapa obsahovat
 
dané kategorie. Ve&nbsp;svém důsledku bude reklasifikovaná, výstupní rastrová mapa obsahovat
celkem tři kategorie (tj.&nbsp;hodnoty). Reklasifikaci provedeme příkazem, předpokládejme, že
+
celkem tři kategorie (tj.&nbsp;hodnoty). Reklasifikaci provedeme příkazem (předpokládejme, že
máme reklasifikační tabulku předpřipravenu v&nbsp;textovém souboru:
+
máme reklasifikační tabulku předpřipravenu v&nbsp;textovém souboru):
  
<pre>
+
#reklasifikace rastrových dat
#reklasifikace rastrových dat
+
#
#
+
GRASS:~&nbsp;>&nbsp;r.reclass input=temp1 output=r_ndvi < reklas-ndvi.txt
GRASS:~&nbsp;>&nbsp;r.reclass input=temp1 output=r_ndvi < reklas-ndvi.txt</tr>
 
</pre>
 
  
 
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://grass.fsv.cvut.cz/~martin/vyuka/zzoz/cv4/tb-r_ndvi.txt tb-r_ndvi.txt]):
+
([http://gama.fsv.cvut.cz/~landa/geowikicz_data/zod/cv4/tb-r_ndvi.txt tb-r_ndvi.txt]):
  
<pre>
+
1 red
1 red
+
2 yellow
2 yellow
+
3 0 136 26
3 0 136 26
 
</pre>
 
  
<pre>
+
#nastavení vhodné tabulky barev
#nastavení vhodné tabulky barev
+
#
#
+
GRASS:~&nbsp;>&nbsp;r.colors r_ndvi col=rules < tb-r_ndvi.txt
GRASS:~&nbsp;>&nbsp;r.colors r_ndvi col=rules < tb-r_ndvi.txt
+
#
#zobrazení reklasifikovaných dat
+
#zobrazení reklasifikovaných dat
GRASS:~&nbsp;>&nbsp;d.rast r_ndvi
+
#
GRASS:~&nbsp;>&nbsp;d.rast.leg r_ndvi
+
GRASS:~&nbsp;>&nbsp;d.rast r_ndvi
</pre>
+
GRASS:~&nbsp;>&nbsp;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 645: Řádek 626:
 
A ještě krátkou statistiku reklasifikované mapy:
 
A ještě krátkou statistiku reklasifikované mapy:
  
<pre>
+
GRASS:~&nbsp;>&nbsp;r.report map=r_ndvi units=c,p,h
GRASS:~&nbsp;>&nbsp;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|
+-----------------------------------------------------------------------------+
+
+-----------------------------------------------------------------------------+
</pre>
 
  
 
== praktická úloha: vytvoření mapy s teplotou povrchu ==
 
== praktická úloha: vytvoření mapy s teplotou povrchu ==
Řádek 679: Řádek 659:
 
v&nbsp;aktuálním mapsetu (student) pod jménem <u>tm6_null</u>:
 
v&nbsp;aktuálním mapsetu (student) pod jménem <u>tm6_null</u>:
  
<pre>
+
#kopie tm6@PERMANENT -> tm6_null@student
#kopie tm6@PERMANENT -> tm6_null@student
+
#
#
+
GRASS:~&nbsp;>&nbsp;g.copy rast=tm6,tm6_null
GRASS:~&nbsp;>&nbsp;g.copy rast=tm6,tm6_null
 
</pre>
 
  
 
Nyní nahradíme hodnotu "0" hodnotou NULL (žádná data) - nulové hodnoty
 
Nyní nahradíme hodnotu "0" hodnotou NULL (žádná data) - nulové hodnoty
Řádek 690: Řádek 668:
 
příkazem:
 
příkazem:
  
<pre>
+
#pro jistotu nastavit aktivní region
#pro jistotu nastavit aktivní region
+
#
#
+
GRASS:~&nbsp;>&nbsp;g.region rast=tm6
GRASS:~&nbsp;>&nbsp;g.region rast=tm6
+
#
#nahrazení hodnoty "0" hodnotou NULL
+
#nahrazení hodnoty "0" hodnotou NULL
#
+
#
GRASS:~&nbsp;>&nbsp;r.null map=tm6_null setnull=0
+
GRASS:~&nbsp;>&nbsp;r.null map=tm6_null setnull=0
</pre>
 
  
 
Dle zadaných vzorců (převzato
 
Dle zadaných vzorců (převzato
 
z&nbsp;[http://mpa.itc.it/grasstutor/index.phtml Open source GIS: A&nbsp;GRASS GIS Approach])
 
z&nbsp;[http://mpa.itc.it/grasstutor/index.phtml Open source GIS: A&nbsp;GRASS GIS Approach])
provedeme [http://grass.fsv.cvut.cz/~martin/vyuka/zzoz/cv4/temp-rovnice.txt sérii] výpočtů:
+
provedeme [http://gama.fsv.cvut.cz/~landa/geowikicz_data/zod/cv4/temp-rovnice.txt sérii] výpočtů:
  
<pre>
+
GRASS:~&nbsp;>&nbsp;r.mapcalc 'tm6rad = 1.238 + (15.6 - 1.238) * tm6_null / 255'
GRASS:~&nbsp;>&nbsp;r.mapcalc 'tm6rad = 1.238 + (15.6 - 1.238) * tm6_null / 255'</tr>
+
GRASS:~&nbsp;>&nbsp;r.mapcalc 'temp_kelvin = 1260.56 / (log (607.76 / tm6rad + 1))'
GRASS:~&nbsp;>&nbsp;r.mapcalc 'temp_kelvin = 1260.56 / (log (607.76 / tm6rad + 1))'</tr>
+
GRASS:~&nbsp;>&nbsp;r.mapcalc 'temp_celsius = temp_kelvin - 273.15'
GRASS:~&nbsp;>&nbsp;r.mapcalc 'temp_celsius = temp_kelvin - 273.15'</tr>
 
</pre>
 
  
 
Před vizualizací výsledné mapy nastavíme vhodnou tabulku barev
 
Před vizualizací výsledné mapy nastavíme vhodnou tabulku barev
([http://grass.fsv.cvut.cz/~martin/vyuka/zzoz/cv4/tb-teplota.txt tb-teplota.txt]):
+
([http://gama.fsv.cvut.cz/~landa/geowikicz_data/zod/cv4/tb-teplota.txt tb-teplota.txt]):
  
<pre>
+
GRASS:~&nbsp;>&nbsp;r.colors map=temp_celsius col=rules < tb-teplota.txt
GRASS:~&nbsp;>&nbsp;r.colors map=temp_celsius col=rules < tb-teplota.txt
+
GRASS:~&nbsp;>&nbsp;d.rast temp_celsius    
GRASS:~&nbsp;>&nbsp;d.rast temp_celsius    
 
</pre>
 
  
 
[[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 723: Řádek 696:
 
budou zajímat pouze parametry <tt>rast</tt> pro rastrová a&nbsp;<tt>vect</tt> pro vektorová data.
 
budou zajímat pouze parametry <tt>rast</tt> pro rastrová a&nbsp;<tt>vect</tt> pro vektorová data.
 
Nic nám tedy nebrání odstranit již nepotřebné rastrové vrstvy, které
 
Nic nám tedy nebrání odstranit již nepotřebné rastrové vrstvy, které
vznikly během výpočtu <font class="mapa">temp_celsius</tt>:
+
vznikly během výpočtu <u>temp_celsius</u>:
 +
 
 +
#odstranění zvolených rastrových mapových vrstev
 +
#
 +
GRASS:~&nbsp;>&nbsp;g.remove rast=tm6_null,tm6rad,temp_kelvin
 +
 
  
<pre>
+
[ [[Zpracování obrazových dat]] ]
#odstranění zvolených rastrových mapových vrstev
 
#
 
GRASS:~&nbsp;>&nbsp;g.remove rast=tm6_null,tm6rad,temp_kelvin
 
</pre>
 

Aktuální verze z 19. 11. 2008, 20:54

[ Zpracování obrazových dat ]

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ů

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.

Tab. č.1: Operátory
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ý -
Tab č.2: Funkce
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.
Tab č.3: Speciální proměnné
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".

Obr č.1 Vstupní kanály tm1 a tm2
Obr č.2: (a) tm1 + tm2; (b) tm1 - tm2; (c) tm1 * tm2; (d) tm1 / tm2; (e) 1.0 * tm1 / tm2

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
Obr č.3: Podíl obrazu (tm7/tm4) - tabulka barev rainbow

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.

Obr č.4: NDVI - tabulka barev rainbow

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
Obr č.5: NDVI - vlastní tabulka barev

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
Obr č.6: Zobrazení NDVI a legendy

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
Obr č.7: Zobrazení NDVI a podruhé

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
Obr č.8: Reklasifikovaný NDVI
Obr č.9: Reklasifikovaný NDVI s legendou

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	  
Obr č.10: Teplota povrchu

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


[ Zpracování obrazových dat ]