153ZODH / 13. cvičení: Porovnání verzí

Z GeoWikiCZ
Skočit na navigaci Skočit na vyhledávání
m (155zddp)
 
(Není zobrazeno 37 mezilehlých verzí od stejného uživatele.)
Řádek 1: Řádek 1:
{{Upravit}}
+
{{Zastaralé|155ZDDP}}
< [[153YZOD Zpracování obrazových dat|Stránky předmětu]] {{bullet}} [[153YZOD Zpracování obrazových dat - cvičení 12|Předchozí cvičení]] {{bullet}} [[153YZOD Zpracování obrazových dat - cvičení 14|Další cvičení]]
+
{{Cvičení|153ZODH|13|Úvod do GRASS/R}}
 
 
<div style="font-size: 13pt; font-weight: bold">Úvod do GRASS/R</div>
 
__TOC__
 
 
== Osnova ==
 
== Osnova ==
  
Řádek 9: Řádek 6:
  
 
=== Seznam příkazů ===
 
=== Seznam příkazů ===
 +
 +
* {{GrassPrikaz|g.mlist}}
 +
* {{GrassPrikaz|i.group}}
 +
* {{GrassPrikaz|g.region}}
 +
* {{GrassPrikaz|r.out.gdal}}
  
 
== Rozhraní GRASS/R ==
 
== Rozhraní GRASS/R ==
  
Základní informace o rozhraní R pro [[GRASS GIS]] jsou dostupné [[GRASS/R|zde]].
+
Základní informace o rozhraní {{freegis|GRASS GIS / Propojení s R|zde}}.
  
 
== Výpočet NDVI ==
 
== Výpočet NDVI ==
  
Spustíme GRASS s lokací '''nc_spm_08''', z příkazové řádky GRASSu spustíme prostředí R a nahrajeme balíček 'spgrass6'.
+
Spustíme GRASS s lokací '''nc_spm_08'''. K dispozici máme snímky {{wikipedia|Landsat 7|lang=en}} z roku 2002.
 
 
GRASS (nc_spm_08):~ > R
 
> library(spgrass6)
 
  
K dispozici máme snímky [http://en.wikipedia.org/wiki/Landsat_7 LandSat 7] z roku 2002.
+
<source lang="bash">
 
+
g.mlist type=rast pattern='lsat7*'
> execGRASS('g.mlist', parameters=list(type='rast', pattern='lsat7*'))
+
</source>
  
 
<pre>
 
<pre>
Řádek 37: Řádek 36:
 
</pre>
 
</pre>
  
 +
Vytvoříme obrazovou skupinu, která obsahuje 1-3, 4-5 a 7 pásmo.
 +
 +
<source lang="bash">
 +
i.group group=lsat7_2002 input=lsat7_2002_10,lsat7_2002_20,lsat7_2002_30,lsat7_2002_40,lsat7_2002_50,lsat7_2002_70
 +
</source>
 +
 +
Nastavíme region a data vyexportujeme do formátu {{wikipedia|GeoTIFF|lang=en}}.
 +
 +
<source lang="bash">
 +
g.region rast=lsat7_2002_10
 +
r.out.gdal input=lsat7_2002 format=GTiff type=Byte output=lsat7_2002.tif
 +
</source>
 +
 +
Soubor by měl obsahovat 6 kanálů.
 +
 +
<source lang="bash">
 +
gdalinfo lsat7_2002.tif | grep Band
 +
</source>
 +
 +
<pre>
 +
Band 1 Block=527x2 Type=Byte, ColorInterp=Gray
 +
Band 2 Block=527x2 Type=Byte, ColorInterp=Undefined
 +
Band 3 Block=527x2 Type=Byte, ColorInterp=Undefined
 +
Band 4 Block=527x2 Type=Byte, ColorInterp=Undefined
 +
Band 5 Block=527x2 Type=Byte, ColorInterp=Undefined
 +
Band 6 Block=527x2 Type=Byte, ColorInterp=Undefined
 +
</pre>
 +
 +
Z příkazové řádky GRASSu spustíme prostředí R a nahrajeme balíček 'spgrass6' a 'RColorBrewer'.
 +
 +
GRASS (nc_spm_08):~ > R
 +
> library(rgdal)
 +
> library(spgrass6)
 +
> library(RColorBrewer)
 +
 +
<source lang="bash">
 +
> GDALinfo('lsat7_2002.tif')
 +
</source>
 +
 +
<pre>
 +
rows        475
 +
columns    527
 +
bands      6
 +
origin.x        629992.5
 +
origin.y        214975.5
 +
res.x      28.5
 +
res.y      28.5
 +
oblique.x  0
 +
oblique.y  0
 +
driver      GTiff
 +
projection  +proj=lcc +lat_1=36.16666666666666 +lat_2=34.33333333333334
 +
+lat_0=33.75 +lon_0=-79 +x_0=609601.22 +y_0=0 +ellps=GRS80 +datum=NAD83
 +
+units=m +no_defs
 +
file        lsat7_2002.tif
 +
apparent band summary:
 +
  GDType Bmin Bmax
 +
1  Byte    0  255
 +
2  Byte    0  255
 +
3  Byte    0  255
 +
4  Byte    0  255
 +
5  Byte    0  255
 +
6  Byte    0  255
 +
</pre>
 +
 +
Vyexportovaná data načteme.
 +
 +
<source lang="bash">
 +
> lsat <- readGDAL('lsat7_2002.tif')
 +
> summary(lsat)
 +
</source>
 +
 +
<pre>
 +
Object of class SpatialGridDataFrame
 +
Coordinates:
 +
      min    max
 +
x 629992.5 645012
 +
y 214975.5 228513
 +
Is projected: TRUE
 +
proj4string :
 +
[+proj=lcc +lat_1=36.16666666666666 +lat_2=34.33333333333334
 +
+lat_0=33.75 +lon_0=-79 +x_0=609601.22 +y_0=0 +ellps=GRS80 +datum=NAD83
 +
+units=m +no_defs +towgs84=0,0,0]
 +
Number of points: 2
 +
Grid attributes:
 +
  cellcentre.offset cellsize cells.dim
 +
x          630006.8    28.5      527
 +
y          214989.8    28.5      475
 +
Data attributes:
 +
    band1            band2            band3          band4     
 +
Min.  : 42.00  Min.  : 28.00  Min.  :  1.0  Min.  :  1.00 
 +
1st Qu.: 72.00  1st Qu.: 55.00  1st Qu.: 41.0  1st Qu.: 81.00 
 +
Median : 80.00  Median : 64.00  Median : 57.0  Median : 93.00 
 +
Mean  : 85.47  Mean  : 70.92  Mean  : 66.7  Mean  : 93.16 
 +
3rd Qu.: 92.00  3rd Qu.: 79.00  3rd Qu.: 81.0  3rd Qu.:106.00 
 +
Max.  :255.00  Max.  :255.00  Max.  :255.0  Max.  :255.00 
 +
    band5          band6     
 +
Min.  :  1.0  Min.  :  1.00 
 +
1st Qu.: 71.0  1st Qu.: 34.00 
 +
Median : 88.0  Median : 52.00 
 +
Mean  : 97.6  Mean  : 61.03 
 +
3rd Qu.:116.0  3rd Qu.: 79.00 
 +
Max.  :255.0  Max.  :255.00 
 +
</pre>
 +
 +
Nahradíme názvy jednotlivých kanálů.
 +
 +
<source lang="bash">
 +
> names(lsat)
 +
</source>
 +
 +
<pre>
 +
[1] "band1" "band2" "band3" "band4" "band5" "band6"
 +
</pre>
 +
 +
<source lang="bash">
 +
> names(lsat) <- c('blue1', 'green2', 'red3', 'nir4', 'mir5', 'mir7')
 +
</source>
 +
 +
Vypočteme NDVI, viz {{153YZODCv|4|Podíl obrazu|4. cvičení}}.
 +
 +
<source lang="bash">
 +
> lsat$ndvi <- (lsat$nir4 - lsat$red3) / (lsat$nir4 + lsat$red3)
 +
> summary(lsat$ndvi)
 +
</source>
 +
 +
<pre>
 +
  Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
 +
-0.9565  0.0000  0.2152  0.1900  0.4140  0.9787
 +
</pre>
 +
 +
Nastavíme vhodnou tabulku barev.
 +
 +
<source lang="bash">
 +
> mypal <- brewer.pal(5, "Greens")
 +
> greens <- colorRampPalette(mypal)
 +
</source>
 +
 +
Data zobrazíme.
 +
 +
<source lang="bash">
 +
> image(lsat, "ndvi", col = greens(20))
 +
</source>
  
----
+
{{fig|r-ndvi|NDVI (GRASS/R)|size=640}}
< [[153YZOD Zpracování obrazových dat|Stránky předmětu]] {{bullet}} [[153YZOD Zpracování obrazových dat - cvičení 12|Předchozí cvičení]] {{bullet}} [[153YZOD Zpracování obrazových dat - cvičení 14|Další cvičení]]
 
  
 
{{ZOD}}
 
{{ZOD}}
 +
{{GRASS}}

Aktuální verze z 5. 9. 2014, 19:49

předchozí cvičenístránky předmětudalší cvičení

Úvod do GRASS/R

Osnova

Toto cvičení je zaměřeno propojení GRASS GIS a prostředí pro statistickou analýzu dat R.

Seznam příkazů

Rozhraní GRASS/R

Základní informace o rozhraní zde.

Výpočet NDVI

Spustíme GRASS s lokací nc_spm_08. K dispozici máme snímky Landsat 7 z roku 2002.

g.mlist type=rast pattern='lsat7*'
lsat7_2002_10
lsat7_2002_20
lsat7_2002_30
lsat7_2002_40
lsat7_2002_50
lsat7_2002_61
lsat7_2002_62
lsat7_2002_70
lsat7_2002_80

Vytvoříme obrazovou skupinu, která obsahuje 1-3, 4-5 a 7 pásmo.

i.group group=lsat7_2002 input=lsat7_2002_10,lsat7_2002_20,lsat7_2002_30,lsat7_2002_40,lsat7_2002_50,lsat7_2002_70

Nastavíme region a data vyexportujeme do formátu GeoTIFF.

g.region rast=lsat7_2002_10
r.out.gdal input=lsat7_2002 format=GTiff type=Byte output=lsat7_2002.tif

Soubor by měl obsahovat 6 kanálů.

gdalinfo lsat7_2002.tif | grep Band
Band 1 Block=527x2 Type=Byte, ColorInterp=Gray
Band 2 Block=527x2 Type=Byte, ColorInterp=Undefined
Band 3 Block=527x2 Type=Byte, ColorInterp=Undefined
Band 4 Block=527x2 Type=Byte, ColorInterp=Undefined
Band 5 Block=527x2 Type=Byte, ColorInterp=Undefined
Band 6 Block=527x2 Type=Byte, ColorInterp=Undefined

Z příkazové řádky GRASSu spustíme prostředí R a nahrajeme balíček 'spgrass6' a 'RColorBrewer'.

GRASS (nc_spm_08):~ > R
> library(rgdal)
> library(spgrass6)
> library(RColorBrewer)
> GDALinfo('lsat7_2002.tif')
rows        475 
columns     527 
bands       6 
origin.x        629992.5 
origin.y        214975.5 
res.x       28.5 
res.y       28.5 
oblique.x   0 
oblique.y   0 
driver      GTiff 
projection  +proj=lcc +lat_1=36.16666666666666 +lat_2=34.33333333333334
+lat_0=33.75 +lon_0=-79 +x_0=609601.22 +y_0=0 +ellps=GRS80 +datum=NAD83
+units=m +no_defs 
file        lsat7_2002.tif 
apparent band summary:
  GDType Bmin Bmax
1   Byte    0  255
2   Byte    0  255
3   Byte    0  255
4   Byte    0  255
5   Byte    0  255
6   Byte    0  255

Vyexportovaná data načteme.

> lsat <- readGDAL('lsat7_2002.tif')
> summary(lsat)
Object of class SpatialGridDataFrame
Coordinates:
       min    max
x 629992.5 645012
y 214975.5 228513
Is projected: TRUE 
proj4string :
[+proj=lcc +lat_1=36.16666666666666 +lat_2=34.33333333333334
+lat_0=33.75 +lon_0=-79 +x_0=609601.22 +y_0=0 +ellps=GRS80 +datum=NAD83
+units=m +no_defs +towgs84=0,0,0]
Number of points: 2
Grid attributes:
  cellcentre.offset cellsize cells.dim
x          630006.8     28.5       527
y          214989.8     28.5       475
Data attributes:
     band1            band2            band3           band4       
 Min.   : 42.00   Min.   : 28.00   Min.   :  1.0   Min.   :  1.00  
 1st Qu.: 72.00   1st Qu.: 55.00   1st Qu.: 41.0   1st Qu.: 81.00  
 Median : 80.00   Median : 64.00   Median : 57.0   Median : 93.00  
 Mean   : 85.47   Mean   : 70.92   Mean   : 66.7   Mean   : 93.16  
 3rd Qu.: 92.00   3rd Qu.: 79.00   3rd Qu.: 81.0   3rd Qu.:106.00  
 Max.   :255.00   Max.   :255.00   Max.   :255.0   Max.   :255.00  
     band5           band6       
 Min.   :  1.0   Min.   :  1.00  
 1st Qu.: 71.0   1st Qu.: 34.00  
 Median : 88.0   Median : 52.00  
 Mean   : 97.6   Mean   : 61.03  
 3rd Qu.:116.0   3rd Qu.: 79.00  
 Max.   :255.0   Max.   :255.00  

Nahradíme názvy jednotlivých kanálů.

> names(lsat)
[1] "band1" "band2" "band3" "band4" "band5" "band6"
> names(lsat) <- c('blue1', 'green2', 'red3', 'nir4', 'mir5', 'mir7')

Vypočteme NDVI, viz 4. cvičení.

> lsat$ndvi <- (lsat$nir4 - lsat$red3) / (lsat$nir4 + lsat$red3)
> summary(lsat$ndvi)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
-0.9565  0.0000  0.2152  0.1900  0.4140  0.9787

Nastavíme vhodnou tabulku barev.

> mypal <- brewer.pal(5, "Greens")
> greens <- colorRampPalette(mypal)

Data zobrazíme.

> image(lsat, "ndvi", col = greens(20))
NDVI (GRASS/R)