NEON EDUCATION bio photo

NEON EDUCATION

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

View All Tutorials

This tutorial is a part of a series!

Click below to view all lessons in the series!

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

About

In this tutorial, we will walk through how to remove parts of a raster based on pixel values using a mask from an analysis.

A mask raster layer is a layer that contains pixels that won’t be used in the analysis. In R, these pixels as assigned an NA value.

Raster Masks

Read more about raster masks in R.

First, let’s load the required libraries.

# load libraries
library(raster)
library(rhdf5)
library(rgdal)

# be sure to set your working directory
# setwd("~/Documents/data/NEONDI-2016") # Mac
# setwd("~/data/NEONDI-2016")  # Windows

## import functions
# install devtools (only if you have not previously intalled it)
#install.packages("devtools")
# call devtools library
#library(devtools)

# install from github
#install_github("lwasser/neon-aop-package/neonAOP")
# call library
library(neonAOP)


#source("/Users/lwasser/Documents/GitHub/neon-aop-package/neonAOP/R/aop-data.R")

Import LiDAR data

To begin, we will open the NEON LiDAR Digital Surface and Digital Terrain Models (DSM and DTM) which are in Geotiff format.

# import aspect data from previous lesson
teak_nsAspect <- raster("outputs/TEAK/TEAK_nsAspect.tif")

# North facing slope = 1
# South facing slope = 2

# legend outside of the plot region
# make room for a legend
par(xpd = FALSE, mar=c(5.1, 4.1, 4.1, 4.5))

plot(teak_nsAspect, 
     col=c("white","blue","green"),
     main="North and South Facing Slopes \n Lower Teakettle",
     legend=F)

# allow legend to plot outside of bounds
par(xpd=TRUE)

legend((par()$usr[2] + 20), 4103300, # set xy legend location
       legend = c("North", "South"),
       fill = c("blue", "green"), 
       bty="n") # turn off border

Mask Data

Once we have created a threhold classified raster, we can use it for different things. One application is to use it as an analysis mask for another dataset.

Let’s try to find all pixels that have an NDVI value >.6 and are north facing.

# open NEON NDVI data
ndvi <- raster("NEONdata/D17-California/TEAK/2013/spectrometer/veg_index/TEAK_NDVI.tif")
ndvi

## class       : RasterLayer 
## dimensions  : 577, 543, 313311  (nrow, ncol, ncell)
## resolution  : 1, 1  (x, y)
## extent      : 325963, 326506, 4102905, 4103482  (xmin, xmax, ymin, ymax)
## coord. ref. : +proj=utm +zone=11 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0 
## data source : /Users/lwasser/Documents/data/1_NEON-DI2016/NEONdata/D17-California/TEAK/2013/spectrometer/veg_index/TEAK_NDVI.tif 
## names       : TEAK_NDVI 
## values      : -0.2252816, 0.9218391  (min, max)

hist(ndvi,
     main="NDVI for Lower Teakettle Field Site")

## Warning in .hist1(x, maxpixels = maxpixels, main = main, plot = plot, ...):
## 32% of the raster cells were used. 100000 values used.

# let's create a mask
ndvi[ndvi<.6] <- NA
plot(ndvi,
     main="NDVI > .6")

n.face.aspect <- teak_nsAspect==1

# mask out only pixels that are north facing and NDVI >.6
nFacing.ndvi <- mask(n.face.aspect, ndvi)

plot(nFacing.ndvi,
     main="North Facing Locations \n NDVI > .6",
     legend=F)

Export Classified Raster

# export geotiff 
writeRaster(nFacing.ndvi,
            filename="outputs/TEAK/TEAK_n_ndvi6.tif",
            format="GTiff",
            options="COMPRESS=LZW",
            overwrite = TRUE,
            NAflag = -9999)

Get Lesson Code

(some browsers may require you to right click.)