#' Contour plot of pollutant concentration
#'
#' \code{contourPlot} plots a contour map of pollutants. This function has been
#' deprecated since version 2.0.0 and will be removed very soon.
#' Use \code{\link{contourPlot2}}.
#'
#' This is a convenience function to plot contour levels of a pollutant matrix
#' with \code{ggplot2}.
#'
#' @param data A dataframe containing data to be plotted in the form of X, Y and
#' Z (levels).
#' @param domain An array with min X, max X, min Y, max Y, number of ticks on X
#' axis, number of ticks on Y axis (optional).
#' @param background String containing the path to the png file to be plotted as
#' a basemap (optional).
#' @param underlayer Array of strings containing layers to be plotted between
#' basemap and contour plot (optional).
#' @param overlayer Array of strings containing layers to be plotted on top of
#' the contour plot (optional).
#' @param legend (string) Legend title (optional).
#' @param levels Array of levels for contour plot. If not set, automatic levels
#' are plotted.
#' @param transparency float (between 0 and 1, default=0.66). Transparency level
#' of the contour plot.
#' @param smoothness integer factor to improve the horizontal resolution
#' (smaller cells) by bilinear interpolation.
#' @param colors Color palette for contour plot
#' @param bare Boolean (default FALSE) parameter to completely remove axis,
#' legend, titles and any other graphical element from the plot.
#' @param size float with the thickness of the contour line.
#' @param cover boolean (default TRUE) to specify whether the contour plot
#' should be filled or not.
#'
#' @return A \code{ggplot2} plot.
#'
#' @examples
#' \dontrun{
#' # Load example data in long format
#' data(volcano)
#' volcano3d <- reshape2::melt(volcano)
#' names(volcano3d) <- c("x", "y", "z")
#' # Contour plot with default options
#' contourPlot(volcano3d)
#'
#' # Import variable CONCAN from inpufile, convert km to m (k = 1000):
#' data <- importRaster(paste0(dir, inputfile),
#' k = 1000,
#' variable = "CONCAN")
#'
#' # Simple contour plot
#' contourPlot(data)
#'
#' # Specifiy (sub)domain to be plotted; background image; legend title and
#' # pollutant levels.
#' contourPlot(data,
#' domain(500000, 510000, 6000000, 6010000, 7, 7),
#' background = "img/background.png",
#' legend = "no2 [ug/m3]",
#' levels = c(10, 20, 30, 40))
#'
#' # Add underlayer (same for overlayer)
#' library(ggplot2)
#' library(maptools)
#' perimetro <- readShapeLines("path_to/perimetro.shp")
#' perimetro <- fortify(perimetro)
#' strada <- readShapeLines("path_to/strada.shp")
#' strada <- fortify(strada)
#' myUnderlayer <- vector(mode = "list", length = 2)
#' myUnderlayer[[1]] <- geom_polygon(data = perimetro,
#' aes(long, lat, group = group),
#' colour = "black",
#' fill = NA,
#' size = 0.1,
#' alpha = 0.5)
#' myUnderlayer[[2]] <- geom_path(data = strada,
#' aes(long, lat, group = group),
#' colour = "grey",
#' size = 0.1,
#' alpha = 0.5)
#' contourPlot(data = test,
#' background = "path_to/basemap.png",
#' underlayer = myUnderlayer)
#'
#' # Change default colour palette
#' contourPlot(data = test,
#' colors = RColorBrewer::brewer.pal(3, name = "PiYG"))
#' }
#'
#' @import ggplot2
#' @importFrom grDevices colorRampPalette
#'
#' @export
#'
contourPlot <- function(data,
domain = NULL,
background = NULL,
underlayer = NULL,
overlayer = NULL,
legend = NULL,
levels = NULL,
size = 0.,
cover = TRUE,
transparency = 0.66,
smoothness = 1.,
colors = NULL,
bare = FALSE) {
.Deprecated("contourPlot2")
if (!requireNamespace("raster", quietly = TRUE)) {
stop("Please install raster from CRAN.", call. = FALSE)
}
# Local binding for variables
x <- y <- z <- NULL
# Convert input to raster
tt <- raster::rasterFromXYZ(data)
# Resample raster
tt <- raster::disaggregate(tt, fact = smoothness, method = 'bilinear')
# Define plot domain
if (missing(domain)) {
xmin <- raster::xmin(tt) # x coordinates minimum
xmax <- raster::xmax(tt) # x coordinates max
ymin <- raster::ymin(tt) # y coordinates min
ymax <- raster::ymax(tt) # y coordinates max
nx <- 5 # number of ticks along x axis
ny <- 5 # number of ticks along y axis
} else {
xmin <- domain[1] # x coordinates minimum
xmax <- domain[2] # x coordinates max
ymin <- domain[3] # y coordinates min
ymax <- domain[4] # y coordinates max
nx <- domain[5] # number of ticks along x axis
ny <- domain[6] # number of ticks along y axis
}
# Automatic scales
if (missing(levels)) {
nlevels <- 7
levels <- pretty(range(raster::values(tt), na.rm = TRUE),
n = nlevels, min.n = 4)
}
lab_levels <- levels
# Extend data domain to be plotted
for (i in (1:1)) {
for (idx in seq(1, 1)) {
# Extend top and boottom rows
ttx1 <- raster::crop(tt, raster::extent(tt, 1, 1, 1,
raster::ncol(tt)))
raster::ymin(ttx1) <- raster::ymax(tt)
raster::ymax(ttx1) <- raster::ymax(tt) + raster::res(tt)[2]
ttxN <- raster::crop(tt, raster::extent(tt, raster::nrow(tt),
raster::nrow(tt), 1,
raster::ncol(tt)))
raster::ymax(ttxN) <- raster::ymin(tt)
raster::ymin(ttxN) <- raster::ymin(tt) - raster::res(tt)[2]
ttE <- raster::merge(tt, ttx1)
ttE <- raster::merge(ttE, ttxN)
# Extend left and right columns
tty1 <- raster::crop(tt, raster::extent(tt, 1,
raster::nrow(tt), 1, 1))
raster::xmin(tty1) <- raster::xmin(tt) - raster::res(tt)[1]
raster::xmax(tty1) <- raster::xmin(tt)
ttyN <- raster::crop(tt, raster::extent(tt, 1,
raster::nrow(tt),
raster::ncol(tt),
raster::ncol(tt)))
raster::xmin(ttyN) <- raster::xmax(tt)
raster::xmax(ttyN) <- raster::xmax(tt) + raster::res(tt)[1]
ttE <- raster::merge(ttE, tty1)
ttE <- raster::merge(ttE, ttyN)
tt <- ttE
}
et <- raster::extent(raster::xmin(ttE) - 1 * raster::res(tt)[1],
raster::xmax(ttE) + 1 * raster::res(tt)[1],
raster::ymin(ttE) - 1 * raster::res(tt)[2],
raster::ymax(ttE) + 1 * raster::res(tt)[2])
mv <- min(raster::values(tt), na.rm = TRUE)
j <- 1
while (lab_levels[j] < mv) {
j <- j + 1
}
if (j != 1) {
ev <- lab_levels[j - 1]
} else {
ev <- lab_levels[1]
}
# ev <- ev - ev/10.
if (lab_levels[1] < 0.) {
ev <- lab_levels[1] - 1.
}
ttE <- raster::extend(ttE, et, value = ev)
}
# convert raster to dataframe
ttP <- raster::rasterToPoints(ttE)
ttDF <- data.frame(ttP)
colnames(ttDF) <- c("x", "y", "z")
# boundaries of extended domain
xminE <- raster::xmin(ttE)
yminE <- raster::ymin(ttE)
xmaxE <- raster::xmax(ttE)
ymaxE <- raster::ymax(ttE)
# color palette (omit first color)
spectral <- c("#5E4FA2", "#3288BD", "#66C2A5", "#ABDDA4", "#E6F598", "#FFFFBF",
"#FEE08B", "#FDAE61", "#F46D43", "#D53E4F", "#9E0142")
if (is.null(colors)) {
myPalette <- grDevices::colorRampPalette(spectral)
myColors <- myPalette(length(levels) + 1)[-c(1, 1)]
} else {
myPalette <- grDevices::colorRampPalette(colors, alpha = TRUE)
myColors <- myPalette(length(levels))
}
# Legend
if (is.null(legend)) {
legend <- ""
}
# prettify legend title
# TODO: svincolarsi da openair per l'etichetta
if (requireNamespace("openair", quietly = TRUE)) {
lgndname <- openair::quickText(legend, auto.text = TRUE)
} else {
lgndname <- legend
}
# Background image
img <- matrix(data = NA, nrow = 10, ncol = 10)
gimg <- grid::rasterGrob(img, interpolate = FALSE)
if (!missing(background)) {
if (requireNamespace("magick", quietly = TRUE)) {
img <- magick::image_read(background)
gimg <- grid::rasterGrob(img)
} else {
warning("Missing magick package. Please install it to be able to read background basemap.")
}
}
# Underlayer
if (missing(underlayer)) {
underlayer <- geom_blank()
}
# Overlayer
if (missing(overlayer)) {
overlayer <- geom_blank()
}
# If requested, wipe all but main plot
if (isTRUE(bare)) {
baretheme <- theme(axis.line=element_blank(),
axis.text.x=element_blank(),
axis.text.y=element_blank(),
axis.ticks=element_blank(),
axis.title.x=element_blank(),
axis.title.y=element_blank(),
legend.position="none",
panel.background=element_blank(),
panel.border=element_blank(),
panel.grid.major=element_blank(),
panel.grid.minor=element_blank(),
plot.background=element_blank())
} else {
baretheme <- theme()
}
# Contour plot
v <- ggplot(ttDF) +
annotation_custom(gimg, -Inf, Inf, -Inf, Inf) +
underlayer +
stat_hollow_contour(
data = ttDF,
aes(x = x, y = y, z = z, fill = factor(..level..)),
geom = "hollow_polygon",
size = size,
breaks = levels,
alpha = transparency,
cover = cover,
na.rm = TRUE) +
scale_fill_manual(lgndname,
guide = guide_legend(reverse = TRUE, label.vjust = 0),
breaks = levels,
labels = lab_levels,
values = myColors) +
scale_x_continuous(name = "x [m]",
limits = c(xminE, xmaxE),
breaks = seq(xmin, xmax, length.out = nx),
expand = c(0, 0)) +
scale_y_continuous(name = "y [m]",
limits = c(yminE, ymaxE),
breaks = seq(ymin, ymax, length.out = ny),
expand = c(0, 0)) +
theme_bw(base_size = 10, base_family = "sans") +
coord_fixed(ratio = 1, xlim = c(xmin, xmax), ylim = c(ymin, ymax)) +
overlayer +
baretheme
return(v)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.