NEON EDUCATION bio photo

NEON EDUCATION

Devoted to open data and open source in science and education.

View All Tutorials

Tags

LiDAR (9)
R programming (70)
Remote Sensing (12)
Data Visualization (4)
Hyperspectral Remote Sensing (7)
Hierarchical Data Formats (HDF5) (24)
Spatial Data & GIS (18)
Time Series (15)
Phenology (7)
Raster Data (8)
Vector Data (6)
Metadata (1)
Git & GitHub (6)
(1) (1) (1)

Tutorial by R Package

dplyr (8)
ggplot2 (17)
h5py (1)
lubridate (time series) (6)
maps (1)
maptools (3)
plyr (2)
raster (32)
rasterVis (raster time series) (3)
rgdal (GIS) (23)
rgeos (5)
rhdf5 (21)
sp (7)
scales (4)
gridExtra (4)
ggtheme (0)
grid (2)
reshape2 (3)
plotly (6)

View ALL Tutorial Series




Twitter Youtube Github


Blog.Roll

R Bloggers

Overview

About

This activity will walk you through the fundamental principles of working with image raster data in R. R Skill Level: Intermediate

Goals / Objectives

After completing this activity, you will know:
  1. How to import multiple image rasters and create a stack of rasters.
  2. How to plot 3 band, RGB images in R.
  3. How to export single band and multiple band image rasters in R.

Things You'll Need To Complete This Lesson

Tools & Libraries To Install

  • R or R studio to write your code.
  • R packages raster, sp and rgdal

Data to Download

Download the raster data:
Download NEON Teaching Data Subset: Field Site Spatial Data
These remote sensing data files provide information on the vegetation at the National Ecological Observatory Network's San Joaquin Experimental Range and Soaproot Saddle field sites. This data is intended for educational purposes, for access to all the data for research purposes visit the NEON Data Portal.

Recommended Pre-Lesson Reading

Review - Raster Data

Raster or “gridded” data are data that are saved in pixels. In the spatial world, each pixel represents an area on the Earth’s surface. An color image raster, is a bit different from other rasters in that it has multiple bands. Each band represents reflectance values for a particular color or type of light. If the image is RGB, then the bands are in the red, green and blue portions of the spectrum. These colors together create what we know as a full color image.

A color image at the NEON San Joachin site in california. Each pixel in the image represents reflectance in the red, green and blue portions of the electromagnetic spectrum.

Bands & Images in R

Working with multiple rasters using Raster Stacks and Raster Bricks

In a previous post, we loaded a single raster into R. We’ve also made sure we knew the CRS (coordinate reference system) and extent of the dataset among other key metadata attributes. Next, let’s create a raster stack from 3 raster images.

A raster stack is a collection of raster layers. Each raster layer in the stack needs to be in the same projection (CRS), spatial extent and resolution. You might use raster stacks for different reasons. For instance, you might want to group a time series of rasters representing precipitation or temperature into one R object. Or, you might want to create a color images from red, green and blue band derived rasters.

In this lesson, we will stack 3 bands from a multi-band image together to create a final RGB image.

First let’s load up the libraries that we need: sp and raster. To install the raster library you can use install.packages(‘raster’). When you install the raster library, sp should also install. Also install the rgdal library” install.packages(‘rgdal’).

#load the raster, sp, and rgdal packages
library(raster)
library(sp)
library(rgdal)
#Set your working directory to the folder where your data for this workshop
#are stored. NOTE: if you created a project file in R studio, then you don't
#need to set the working directory as it's part of the project.
setwd("~/yourWorkingDirectoryHere")  

Next, let’s create a raster stack. The difficult way to do this is to load our rasters one at a time. But that takes a good bit of effort.

#import tiffs
band19 <- "rasterLayers_tif/band19.tif"
band34 <- "rasterLayers_tif/band34.tif"
band58 <- "rasterLayers_tif/band58.tif"

We can also use the list.files command to grab all of the files in a directory. We can use full.names=TRUE to ensure that R will store the directory path in our list of rasters.

#create list of files to make raster stack
rasterlist <-  list.files('rasterLayers_tif', full.names=TRUE)

NOTE: If your list of rasters is located in your main R working directory, you can achieve the same results as above by looking for all files with a ‘.tif’ extension: rasterlist <- list.files('rasterLayers_tif', full.names=TRUE, pattern="tif").

#create raster stack
rgbRaster <- stack(rasterlist)

#check to see that you've created a raster stack and plot the layers
rgbRaster
plot(rgbRaster)
All rasters in the rasterstack plotted.

Plot an RGB Image

You can plot an RGB image from a raster stack. You need to specify the order of the bands when you do this. In our raster stack, band 19 which is the blue band, is first in the list, whereas band 58 which is the red band is last. Thus the order is 3,2,1 to ensure the red band is rendered first as red.

#plot an RGB version of the stack
plotRGB(rgbRaster,r=3,g=2,b=1, scale=800, stretch = "Lin")
To plot an RGB image in R, you need to specify which rasters to render on the red, green and blue bands.

Explore Raster Values - Histograms

You can also explore the data.

#look at histogram of reflectance values for all rasters
hist(rgbRaster)
Histogram of reflectance values for each raster in the raster stack.
#remember that crop function? You can crop all rasters within a raster stack too
#finally you can crop all rasters within a raster stack!
rgbCrop <- c(256770.7,256959,4112140,4112284)
rgbRaster_crop <- crop(rgbRaster, rgbCrop)
plot(rgbRaster_crop)
Cropped rasters in the raster stack.

Raster Bricks in R

Now we have a list of rasters in a stack. These rasters are all the same extent CRS and resolution but a raster brick will create one raster object in R that contains all of the rasters we can use this object to quickly create RGB images. Raster bricks are more efficient objects to use when processing larger datasets. This is because the computer doesn’t have to spend energy finding the data - it is contained within the object.

#create raster brick
RGBbrick <- brick(rgbRaster)
plotRGB(RGBbrick,r=3,g=2,b=1, scale=800, stretch = "Lin")

While the brick might seem similar to the stack, we can see that it’s very different when we look at the size of the object.

#the brick contains ALL of the data stored in one object
#the stack contains links or references to the files stored on your computer
object.size(rgbBrick)

Output

> 5759448 bytes	

View the size of the stack:

object.size(rgbRaster)

Output:

> 40408 bytes

Write a raster to a Geotiff File in R

We can write out the raster in tiff format as well. When we do this it will copy the CRS, extent and resolution information so the data will read properly into a GIS as well. Note that this writes the raster in the order they are in - in the stack. In this case, the blue (band 19) is first but it’s looking for the red band first (RGB). One way around this is to generate a new raster stack with the rasters in the proper order - red, green and blue format.

#Make a new stack in the order we want the data in: 
finalRGBstack <- stack(rgbRaster$band58,rgbRaster$band34,rgbRaster$band19)
#write the geotiff - change overwrite=TRUE to overwrite=FALSE if you want to make sure you don't overwrite your files!
writeRaster(finalRGBstack,"rgbRaster.tif","GTiff", overwrite=TRUE)

Import A Multi-Band Image into R

You can import a multi-band image into R too. To do this, you import the file as a stack rather than a raster (which brings in just one band). Let’s import the raster than we just created above.

#Import Multi-Band raster
newRaster <- stack("rgbRaster.tif") 
#then plot it
plot(newRaster)
plotRGB(newRaster,r=1,g=2,b=3, scale=800, stretch="lin")

Challenge

  1. A color infrared image is a combination of the gree, red and near-infrared bands. In our case gree <- band 34, red <- band 58 and near-infrared <- band90. Using the same band90 raster that you fixed in step one above, create a color infrared image. HINT: when you use the plotRGB function, you’ll map the bands as follows: blue=band34, green=band58, red=band90. Export your results as a geotiff.

Get Lesson Code

(some browsers may require you to right click.)