Carga de librerías necesarias:

library(raster)
## Loading required package: sp
library(rgdal)
## rgdal: version: 1.2-16, (SVN revision 701)
##  Geospatial Data Abstraction Library extensions to R successfully loaded
##  Loaded GDAL runtime: GDAL 2.1.2, released 2016/10/24
##  Path to GDAL shared files: /usr/share/gdal/2.1
##  GDAL binary built with GEOS: TRUE 
##  Loaded PROJ.4 runtime: Rel. 4.9.3, 15 August 2016, [PJ_VERSION: 493]
##  Path to PROJ.4 shared files: (autodetected)
##  Linking to sp version: 1.2-7

Consultar información de un archivo raster:

GDALinfo("/home/jorge/Documentos/Maestria/PRA/Dataset_GEE/DEM_Bosavita_10m.tif")
## rows        1808 
## columns     1732 
## bands       1 
## lower left origin.x        -73.56242 
## lower left origin.y        5.256726 
## res.x       1.078094e-05 
## res.y       1.07873e-05 
## ysign       -1 
## oblique.x   0 
## oblique.y   0 
## driver      GTiff 
## projection  +proj=longlat +datum=WGS84 +no_defs 
## file        /home/jorge/Documentos/Maestria/PRA/Dataset_GEE/DEM_Bosavita_10m.tif 
## apparent band summary:
##    GDType hasNoDataValue NoDataValue blockSize1 blockSize2
## 1 Float64           TRUE           0          1       1732
## apparent band statistics:
##       Bmin     Bmax   Bmean      Bsd
## 1 2760.348 3045.331 2914.69 50.75652
## Metadata:
## AREA_OR_POINT=Area

Guardar los datos del archivo raster en R:

Data_Bosavita_DSM <- capture.output(
  GDALinfo("/home/jorge/Documentos/Maestria/PRA/Dataset_GEE/DEM_Bosavita_10m.tif")
)

Abrir un archivo raster en R

DSM_BVTA <- 
  raster("/home/jorge/Documentos/Maestria/PRA/Dataset_GEE/DEM_Bosavita_10m.tif")

DSM_BVTA
## class       : RasterLayer 
## dimensions  : 1808, 1732, 3131456  (nrow, ncol, ncell)
## resolution  : 1.078094e-05, 1.07873e-05  (x, y)
## extent      : -73.56242, -73.54375, 5.256726, 5.27623  (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0 
## data source : /home/jorge/Documentos/Maestria/PRA/Dataset_GEE/DEM_Bosavita_10m.tif 
## names       : DEM_Bosavita_10m 
## values      : 2760.348, 3045.331  (min, max)

Ver estadísticas descriptivas del raster cargado:

summary(DSM_BVTA)
## Warning in .local(object, ...): summary is an estimate based on a sample of 1e+05 cells (3.19% of all cells)
##         DEM_Bosavita_10m
## Min.            2759.373
## 1st Qu.         2893.265
## Median          2911.083
## 3rd Qu.         2941.627
## Max.            3045.274
## NA's               0.000

Para calcular las estadíticas sobre todo el raster y no sobre una muetra aleatoria se puede usar:

summary(DSM_BVTA, maxsamp = ncell(DSM_BVTA))
##         DEM_Bosavita_10m
## Min.            2759.258
## 1st Qu.         2893.138
## Median          2911.093
## 3rd Qu.         2941.427
## Max.            3045.343
## NA's         1545995.000

Para visualizar la capa en ggplot2 se debe convertir el raster en un dataframe:

DSM_BVTA_df <- raster::as.data.frame(DSM_BVTA, xy = TRUE)

Ahora se puede ver la imagen como un dataframe:

str(DSM_BVTA_df)
## 'data.frame':    3131456 obs. of  3 variables:
##  $ x               : num  -73.6 -73.6 -73.6 -73.6 -73.6 ...
##  $ y               : num  5.28 5.28 5.28 5.28 5.28 ...
##  $ DEM_Bosavita_10m: num  NA NA NA NA NA NA NA NA NA NA ...

Visualización de la imagen usango ggplot():

library(ggspatial)
## Loading required package: ggplot2
ggplot() +
  geom_spatial(data = DSM_BVTA, aes(fill = band1)) +
  scale_fill_viridis_c() +
  coord_equal()

Ver el sistema de coordenadas de la imagen:

crs(DSM_BVTA)
## CRS arguments:
##  +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0

Calcular valores mínimos y máximos del raster:

DSM_BVTA <- setMinMax(DSM_BVTA)

Valor mínimo:

minValue(DSM_BVTA)
## [1] 2759.258

Valor máximo:

maxValue(DSM_BVTA)
## [1] 3045.343

Crear un histograma del raster:

library(ggplot2)
ggplot(DSM_BVTA_df, aes(DEM_Bosavita_10m)) +
  geom_histogram()
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## Warning: Removed 1545995 rows containing non-finite values (stat_bin).

#ggplot() + geom_histogram(data = DSM_BVTA_df)

Cambiar el ancho de banda con la que se muestra el histograma:

ggplot(DSM_BVTA_df, aes(DEM_Bosavita_10m)) +
  geom_histogram(binwidth = 1)
## Warning: Removed 1545995 rows containing non-finite values (stat_bin).