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

Z GeoWikiCZ
← Porovnání se starší verzí
VizuálníWikitext
m 155zddp
 
(Není zobrazeno 45 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 10: Řádek 7:
=== Seznam příkazů ===
=== Seznam příkazů ===


== Úvod do R ==
* {{GrassPrikaz|g.mlist}}
* {{GrassPrikaz|i.group}}
* {{GrassPrikaz|g.region}}
* {{GrassPrikaz|r.out.gdal}}
 
== Rozhraní GRASS/R ==
 
Základní informace o rozhraní {{freegis|GRASS GIS / Propojení s R|zde}}.


[[Image:Rlogo.png|100px|left]]
== Výpočet NDVI ==
[http://www.r-project.org R] je programovací jazyk a prostředí určené pro ''statistickou analýzu dat a jejich grafické zobrazení''. Jde o implementaci programovacího jazyka [http://en.wikipedia.org/wiki/S_(programming_language) S] pod svobodnou licencí. R již předstihlo počtem uživatelů komerční S a stalo se faktickým standardem v řadě oblastí statistiky.


Používá se z příkazové řádky, existuje však několik frontendů s grafickým rozhraním jako [http://en.wikipedia.org/wiki/RKWard RKWard], [http://en.wikipedia.org/wiki/R_Commander R Commander] nebo rozšíření do OpenOffice.org Calcu [http://extensions.services.openoffice.org/project/R4Calc R4Calc].
Spustíme GRASS s lokací '''nc_spm_08'''. K dispozici máme snímky {{wikipedia|Landsat 7|lang=en}} z roku 2002.


R bývá také propojováno či využíváno v komerčních softwarech, např. v prostředí [http://cs.wikipedia.org/wiki/PASW PASW] mohou uživatelé přímo psát a spouštět programy v jazyce R nad otevřenými daty.
<source lang="bash">
g.mlist type=rast pattern='lsat7*'
</source>


Další informace na [http://en.wikipedia.org/wiki/R_(programming_language) wikipedii].
<pre>
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
</pre>


=== Balíčky pro analýzu geoprostorových dat ===
Vytvoříme obrazovou skupinu, která obsahuje 1-3, 4-5 a 7 pásmo.


R nabízí až více než [http://cran.mirroring.de/web/packages/ 2000] různých ''rozšíření'' (tzv. balíčku) specializovaných na různé analýzy dat. My se zaměříme na balíčky určené pro práci s geoprostorovými (většinou obrazovými) daty.
<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>


* '''[http://cran.mirroring.de/web/packages/sp/index.html sp]''' - základní balíček definující třídy a metody pro práci s prostorovými daty
Nastavíme region a data vyexportujeme do formátu {{wikipedia|GeoTIFF|lang=en}}.


* [http://cran.mirroring.de/web/packages/maptools/index.html maptools] - manipulace s prostorovými objekty
<source lang="bash">
g.region rast=lsat7_2002_10
r.out.gdal input=lsat7_2002 format=GTiff type=Byte output=lsat7_2002.tif
</source>


* [http://cran.mirroring.de/web/packages/maps/index.html maps] - vykreslování map
Soubor by měl obsahovat 6 kanálů.


* [http://cran.mirroring.de/web/packages/spatial/index.html spatial] - Kriging
<source lang="bash">
gdalinfo lsat7_2002.tif | grep Band
</source>


* [http://cran.mirroring.de/web/packages/spatstat/index.html spatstat] - analýza bodových výskytů
<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>


* [http://cran.mirroring.de/web/packages/splancs/index.html splancs] - analýza bodových výskytů v čase
Z příkazové řádky GRASSu spustíme prostředí R a nahrajeme balíček 'spgrass6' a 'RColorBrewer'.


* [http://cran.mirroring.de/web/packages/spdep/index.html spdep] - autokorelace prostorových objektů
GRASS (nc_spm_08):~ > R
> library(rgdal)
> library(spgrass6)
> library(RColorBrewer)


* [http://cran.mirroring.de/web/packages/gstat/index.html gstat] - geostatistické modelování
<source lang="bash">
> GDALinfo('lsat7_2002.tif')
</source>


* [http://cran.mirroring.de/web/packages/geoR/index.html geoR] - analýza geostatistických dat
<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>


* [http://cran.mirroring.de/web/packages/fields/index.html fields] - analýza geoprostorových dat
Vyexportovaná data načteme.


* [http://cran.mirroring.de/web/packages/spatialCovariance/index.html] - výpočet kovariační matice prostorových dat
<source lang="bash">
> lsat <- readGDAL('lsat7_2002.tif')
> summary(lsat)
</source>


* [http://cran.mirroring.de/web/packages/RArcInfo/index.html RArcInfo] - import dat ve formátu ArcInfo Coverage
<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>


* [http://cran.mirroring.de/web/packages/shapefiles/index.html shapefiles] - čtení a zápis ESRI Shapefile
Nahradíme názvy jednotlivých kanálů.


* [http://cran.mirroring.de/web/packages/RColorBrewer/index.html RColorBrewer] - tabulky barev optimalizované pro tématické mapy
<source lang="bash">
> names(lsat)
</source>


* '''[http://cran.mirroring.de/web/packages/spgrass6/index.html spgrass6]''' - rozhraní pro GRASS GIS 6+
<pre>
[1] "band1" "band2" "band3" "band4" "band5" "band6"
</pre>


* [http://cran.mirroring.de/web/packages/rgdal/index.html rgdal] - rozhraní pro knihovnu [[GDAL/OGR|GDAL]]
<source lang="bash">
> names(lsat) <- c('blue1', 'green2', 'red3', 'nir4', 'mir5', 'mir7')
</source>


<center>
Vypočteme NDVI, viz {{153YZODCv|4|Podíl obrazu|4. cvičení}}.
{|class="border"
|+Tab. č. 1: Datové typy definované balíčkem sp
|-
!Datový typ !! Třída !! Rodičovská třída
|-
| body || SpatialPoints || Spatial
|-
| pixely || SpatialPixels || SpatialPoints
|-
| mřížka || SpatialGrid || SpatialPixels
|-
| linie || Line ||
|-
| linie || Lines || Line
|-
| linie || SpatialLines || Spatial, Line
|-
| polygon || Polygon || Line
|-
| polygony || Polygons || Polygon
|-
| polygon || SpatialPolygons || Spatial, Polygon
|-
|}
</center>


== Rozhraní GRASS/R ==
<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>


== Externí odkazy ==
Data zobrazíme.
* [http://www.stanford.edu/~cengel/spatialanthro/archives/R_GRASS_Spatial.pdf R, GRASS, and Spatial Analysis]


<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, 18:49


Tento text se vztahuje k předmětu, který již není součástí aktuálních studijních plánů.
Stránky navazujícího předmětu 155ZDDP jsou dostupné zde.
 


Ú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)
Citováno z „http://geo102.fsv.cvut.cz:8091/index.php?title=153ZODH_/_13._cvičení&oldid=29414“