Import Raster Data

This Notebook will demonstrate how to import raster GIS data into R.

Import a Landsat 8 Image

Because Landsat images are multi-band we import them using the brick():

library(raster)
yose_l8_rst <- raster::brick("./data/yose_l8_20180822_b2345.tif")
yose_l8_rst

class      : RasterBrick 
dimensions : 2555, 2014, 5145770, 4  (nrow, ncol, ncell, nlayers)
resolution : 30, 30  (x, y)
extent     : 246315, 306735, 4152885, 4229535  (xmin, xmax, ymin, ymax)
crs        : +proj=utm +zone=11 +datum=WGS84 +units=m +no_defs 
source     : yose_l8_20180822_b2345.tif 
names      : yose_l8_20180822_b2345.1, yose_l8_20180822_b2345.2, yose_l8_20180822_b2345.3, yose_l8_20180822_b2345.4 
min values :                     7208,                     6100,                     5357,                     4330 
max values :                    29840,                    30960,                    32595,                    30481 

View a summary of the bands:

summary(yose_l8_rst)

summary is an estimate based on a sample of 1e+05 cells (1.94% of all cells)

        yose_l8_20180822_b2345.1 yose_l8_20180822_b2345.2 yose_l8_20180822_b2345.3
Min.                        7411                     6366                     5609
1st Qu.                     8589                     8072                     7749
Median                      9254                     8908                     8965
3rd Qu.                    10419                    10456                    10960
Max.                       26650                    27824                    29578
NA's                           0                        0                        0
        yose_l8_20180822_b2345.4
Min.                        4933
1st Qu.                    12439
Median                     13980
3rd Qu.                    15809
Max.                       27484
NA's                           0

We can tell from the values that Landsat pixel values are 16-bit (0..65536). Remind ourselves of the band combinations for Landsat 8:

Plot it with plotRGB

Plot a Pseudo True Color:

plotRGB(yose_l8_rst, r=3, g=2, b=1)

Error in grDevices::rgb(RGB[, 1], RGB[, 2], RGB[, 3], alpha = alpha, maxColorValue = scale) : 
  color intensity 9182, not in 0:255

Plot a False Color Composite:

plotRGB(yose_l8_rst, r=4, g=3, b=2)

Error in grDevices::rgb(RGB[, 1], RGB[, 2], RGB[, 3], alpha = alpha, maxColorValue = scale) : 
  color intensity 8748, not in 0:255

Plot it with tmap

To overlay the park boundary, we’ll use tmap:

library(sf)
library(tmap)
epsg_utm11n_wgs84 <- 32611
yose_bnd_utm <- sf::st_read(dsn="./data", layer="yose_boundary") %>% 
  st_transform(epsg_utm11n_wgs84)

Reading layer `yose_boundary' from data source 
  `/home/lindangulopez/Desktop/SpatialRWorkshop/data' using driver `ESRI Shapefile'
Simple feature collection with 1 feature and 11 fields
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: -119.8864 ymin: 37.4947 xmax: -119.1964 ymax: 38.18515
Geodetic CRS:  North_American_Datum_1983

tm_shape(yose_l8_rst) +
  tm_rgb(r=4, g=3, b=2, interpolate = FALSE, max.value = 2 ^ 16) +
tm_shape(yose_bnd_utm) +
  tm_borders(col="red", lwd=2) +
tm_layout(legend.show = FALSE)

stars object downsampled to 888 by 1126 cells. See tm_shape manual (argument raster.downsample)