# Plotting rasters

HOME The 'raster' package has a plot method for Raster objects that has many options to customise the plot. Unfortunately the documentation is in two places, and some arguments are not documented at all. Here I'll cover the arguments which you can use when plotting a RasterLayer object.

On this page we'll look at a complete list of the arguments to the plot function, some commonly needed tweaks, and some example plots.

## Arguments for plotting rasters

The main help page for the plot function is at help(plot, "raster"). We're concerned here with the S4 method for signature 'Raster,ANY'. A number of arguments are described here, then we have the 'dots' argument, ..., which can take a huge variety of named arguments, including those accepted by image.plot function in the fields package.

 Argument Default Description These are the documented arguments to the plot method: x (missing) A Raster* object; what follows here applies to a RasterLayer object. y (missing) Ignored if x is a RasterLayer maxpixels 5e5 Maximum number of pixels to plot, if necessary the raster will be passed to sampleRegular for subsampling. col (missing) A colour palette to be used; if missing and x has no 'color table',  rev(terrain.colors(255)) is used alpha NULL Transparency, NULL or 1 = not transparent. colNA NA The background colour, used for pixels with value = NA. Default is transparent. add FALSE If TRUE, the plot will be added to an existing plot. ext NULL Defines a part of the raster to be plotted; see ?extent useRaster TRUE If TRUE, the rasterImage function is used for plotting. Otherwise the image function is used. interpolate FALSE If TRUE, the image will be smoothed. addfun NULL Used when multiple layers are plotted. nc, nr (missing) maxnl 16 main (missing) Main plot title. npretty 0 Number of decimals for axis labels when plotting lon/lat. The plot method calls .rasterImagePlot() which is cloned from imagePlot() in the fields package. It has the following new arguments, most of which are described in help(image.plot, "fields"): legend TRUE If FALSE, no legend bar is displayed. horizontal FALSE If TRUE, the legend bar is displayed below the main plot. legend.shrink 0.5 Amount to shrink the length of the legend bar relative to the full height (or width if horizontal=TRUE) of the plot. legend.width 0.6 Width in characters of the legend bar. legend.mar ifelse(horizontal, 3.1, 5.1) Width in characters of legend margin that has the axis. legend.lab NULL Label for the legend bar; default is no label. graphics.reset FALSE If FALSE (default) the plotting region (plt in par) will not be reset and one can add more information onto the image plot (e.g. using functions such as points or lines.) bigplot NULL Plot coordinates for image plot. smallplot NULL Plot coordinates for legend strip. legend.only FALSE If TRUE just add the legend to a plot (see below). breaks NULL Break points in sorted order to indicate the intervals for assigning the colors. Note that if there are nlevel colors there should be (nlevel+1) breakpoints. lab.breaks NULL If breaks are supplied these are text string labels to put at each break value. axis.args NULL A list of additional arguments for the legend axis. legend.args NULL A list of arguments for a complete specification of the legend label. box TRUE Draws a box around the raster image zlim NULL A vector of length 2 setting limits on the values plotted. zlimcol NULL If NULL, values outside zlim are replaced with NA and plotted according to colNA, otherwise the max or min values. fun NULL Function to transform values before plotting. asp (missing) Aspect ratio, default is to have equal scales for x and y axis. These arguments from plot.default are also respected: xlim, ylim axes TRUE

## Examples

As examples, we'll use data for Switzerland in the unmarked package and for the Bernese Oberland in the AHMbook package:

library(unmarked)
data(Switzerland)
str(Switzerland)
rche <- rasterFromXYZ(Switzerland[1:3]) # elevation
plot(rche)

library(AHMbook)
data(BerneseOberland)
str(BerneseOberland)
rbo <- rasterFromXYZ(BerneseOberland[1:3])
plot(rbo)

Plot Switzerland's elevation and forest cover data with no axes, and using green for forest cover:

rchf <- rasterFromXYZ(Switzerland[c(1,2,4)]) # forest
par(mfrow=1:2)
plot(rchf, col=rampGreens(), axes = FALSE, main="Forest cover")
plot(rche, axes = FALSE, main="Elevation")

There's rather a lot of white space, but not quite enough for the legend on the right plot. We'll reduce the margins and also manually reduce the height of the plotting window.

par(mfrow=1:2, mar=c(0,0,3,1)+0.1)
plot(rchf, col=rampGreens(), axes = FALSE, main="Forest cover")
plot(rche, axes = FALSE, main="Elevation")

Note that the mar settings only affect the image plot, not the legend bar, and that the last value - the right margin width - determines the gap between the image and the legend bar. The figure is more compact but too little space for the legend labels. We can fix that with the legend.mar argument, making it just over 8 characters wide, and we'll remove the boxes too:

plot(rchf, col=rampGreens(), axes = FALSE, main="Forest cover", legend.mar=8, box=FALSE)
plot(rche, axes = FALSE, main="Elevation", legend.mar=8, box=FALSE)

With these plots, I keep the box around the plot, manually resize the plotting window and check replotting, then remove the box for the final version.

Horizontal legend bar

Unless you have axes and axis labels in the bottom margin, you will probably want to reduce the bottom margin, as that is between the main plot and the legend bar.

par(mar=c(0,0,0,0)+0.1)
plot(rchf, col=rampGreens(256), axes = FALSE, horizontal = TRUE, box=FALSE)

Zooming in

We can "zoom" to a smaller area of the raster with either the ext argument (which needs an extent object) or xlim and ylim. I'll zoom the Swiss forest map to the same extent as the Bernese Oberland map on the left and on the right, zoom to the western tip of Switzerland, using black to mask the land outside the border.

par(mfrow=1:2, mar=c(0,0,3,0)+0.1)
plot(rchf, axes = FALSE, main="With ext", ext=extent(rbo), legend=FALSE)
plot(rchf, axes = FALSE, main="With xlim, ylim", colNA = 'black',
xlim=c(910442,960442), ylim=c(53776,103776), legend=FALSE)

Interpolation

These plots look somewhat pixellated. The interpolate argument will smooth the edges.

par(mfrow=1:2, mar=c(0,0,3,0)+0.1)
plot(rbo, axes = FALSE, main="Original", legend=FALSE)
plot(rbo, axes = FALSE, main="Interpolated", interpolate=TRUE, legend=FALSE)

Plots with the same scale and a single legend

Now let's plot elevation for all of Switzerland and the Bernese Oberland together:

par(mfrow=1:2, mar=c(0,0,3,1)+0.1)
plot(rche, axes = FALSE, main="Switzerland")
plot(rbo, axes = FALSE, main="Bernese Oberland")

Quite a bit to do to sort this out! First though I want to have the same colour coding for both plots, even though the second one has a smaller range of elevation values. We can check the min/max values for each raster with:

rche ; rbo

Switzerland's elevation ranges from 190 to 4472m but the Bernese Oberland only from 509 to 3084m. We can use the zlim argument to set the same range for each plot. There should be no values outside the zlim range, so no need to specify zlimcol.

zlim <- range(values(rche), na.rm=TRUE)
plot(rche, axes = FALSE, main="Switzerland", zlim=zlim)
plot(rbo, axes = FALSE, main="Bernese Oberland", zlim=zlim)

Since both plots use the same scale, we only need one legend bar, which we'll  put on the right, making enough space for it with legend.mar. That should be done now, so we remove the boxes.

plot(rche, axes = FALSE, main="Switzerland", zlim=zlim, legend=FALSE, box=FALSE)
plot(rbo, axes = FALSE, main="Bernese Oberland", zlim=zlim, legend.mar=8, box=FALSE)

An issue with having one plot with a legend and one without, is that the plot with the legend is slightly smaller. That doesn't show up above, but look at this:

zlim <- range(values(rbo))
par(mfrow=1:2, mar=c(0,0,3,0.5)+0.1)
plot(rbo, axes = FALSE, main="Before", zlim=zlim, legend=FALSE)
plot(rbo, axes = FALSE, main="After", zlim=zlim, legend.mar=8)

To keep the images the same size, we plot both without the legend, then add the legend later:

plot(rbo, axes = FALSE, main="Before", zlim=zlim, legend=FALSE)
plot(rbo, axes = FALSE, main="After", zlim=zlim, legend=FALSE)
plot(rbo, zlim=zlim, legend.only=TRUE)

Sizes are now ok, but not enough room for the legend labels. Adding legend.mar = 8 makes a mess (try it!) so we'll make space by setting an outer margin:

par(mfrow=1:2, mar=c(0,0,3,0.5)+0.1, oma=c(0.5,0.5,0,2))
plot(rbo, axes = FALSE, main="Before", zlim=zlim, legend=FALSE)
plot(rbo, axes = FALSE, main="After", zlim=zlim, legend=FALSE)
plot(rbo, zlim=zlim, legend.only=TRUE)

Updated 14 July 2019 by Mike Meredith