--- title: 'Operanciones contextuales' author: "Paqui y Pala" date: "Mayo 2017" output: html_document: fig_height: 3 fig_width: 3 number_sections: yes toc: yes pdf_document: fig_height: 3 fig_width: 3 number_sections: yes toc: yes subtitle: Creación de datos espaciales institute: Escuela Internacional de Doctorado (EIDUM) --- # Intro ```{r} library( raster ) ``` # Cambiando la geometría # Pendientes y orientaciones ```{r} mde <- raster( "capasSIG/chuecosMde.tif" ) plot( mde ) pendiente <- terrain( mde ) orientacion <- terrain( mde, "aspect" ) plot( stack( pendiente, orientacion ) ) sombreado <- hillShade( pendiente, orientacion ) plot( sombreado, col = grey.colors( 5 ) ) ``` # Ejemplos de transformación ```{r} r <- raster( ncol = 150, nrow = 150, xmn = -100, xmx = 100, ymn = -100, ymx = 100 ) crs( r ) <- "+proj=lcc +lat_1=48 +lat_2=33 +lon_0=-100 +ellps=WGS84" xy = c(0,0) dist <- distanceFromPoints( r, xy ) plot( dist ) contour( dist / 1000, asp = 1, add = TRUE, nlevels = 10 ) persp( dist, theta = 30, phi = 35, shade = 0.25, col = "green", maxpixels = 1000, expand = 0.5 ) distMax <- 100 - dist persp( distMax, theta = 30, phi = 35, shadow = 15, col = "green", maxpixels = 1000, expand = 0.5 ) distMax[ dist > 50 ] <- 50 plot( distMax ) contour( distMax , asp = 1, add = TRUE, nlevels = 10 ) persp( distMax, theta = 30, phi = 35, shadow = 15, col = "green", maxpixels = 1000, expand = 0.5 ) distMaxSlope <- terrain( distMax ) distMaxAspect <- terrain( distMax, opt = "aspect" ) hsDistMax <- hillShade( slope = distMaxSlope, aspect = distMaxAspect ) plot( hsDistMax, col = grey( seq( 0, 1, 0.01 ) ) ) dist <- distanceFromPoints( r, xy ) contour( dist, asp = 1 ) max <- 50 distMax2 <- sin( acos( dist / max ) ) * max distMax2[ is.nan( distMax2@data@values ) | is.na( distMax2@data@values ) ] <- 0 persp( distMax2, expand = 0.25, maxpixels = 2000, phi = 25, col = "green", border = "darkgreen" ) plot( distMax2 ) contour( distMax2, add = TRUE, nlevels = 25 ) distMaxSlope <- terrain( distMax2 ) distMaxAspect <- terrain( distMax2, opt = "aspect" ) hsDistMax <- hillShade( slope = distMaxSlope, aspect = distMaxAspect ) plot( hsDistMax, col = grey( seq( 0, 1, 0.01 ) ) ) ``` # ¿Qué hay aquí? ```{r, eval=FALSE} click() locator() zoom() ``` ## Para muestrear con puntos ```{r} x <- runif( 100, mde@extent@xmin, mde@extent@xmax ) y <- runif( 100, mde@extent@ymin, mde@extent@ymax ) umCentro <- data.frame( x, y ) plot( mde ) points( umCentro, pch = 4, col = "grey30" ) ``` ```{r} extract( mde, umCentro ) umCentro$cota <- extract( mde, umCentro[ , c( "x", "y" ) ] ) umCentro #umCentro$geologia <- extract( geologia, umCentro ) #umCentro #aggregate( umCentro, list( umCentro$geologia ), median ) ``` ## Para muestreo con polígonos ```{r, eval=FALSE} areaTrabajo <- drawPoly() extract( mde, areaTrabajo, mean ) ``` # Resumen de la sesión ```{r} sessionInfo() ``` # Enlaces * [GIS and Mapping with R](http://www.richardcondit.org/workshops/Rgis/tupper2015/index.html) * [The Visual Raster Cheat Sheet](https://rpubs.com/etiennebr/visualraster) * [GIS and Mapping with R](http://remote-sensing.eu/update-r-spatial-handling-cheatsheet/)