Description Note Note Author(s) References See Also Examples
Functions for visualising 2-D and 3-D data.
Many of the functions are extensions of R's persp or image function.
Other packages that provide visualisation of 3-D data (and which might be better suited) are:
rgl, scatterplot3D, misc3D
.
This package is dedicated to Carlo.
Some of the functions based on persp
will not work properly for all
values of phi
(which turns the plots upside-down).
This is because an assumption is made as to how the perspective plots are viewed.
Karline Soetaert
http://www.rforscience.com/rpackages/visualisation/oceanview/
http://www.rforscience.com/rpackages/visualisation/plot3d/
Functions that are based on the persp function:
persp3D: an extended version of persp.
ribbon3D: a perspective plot as ribbons.
hist3D: 3-D histograms.
scatter3D, points3D, lines3D: colored points, lines, ... in 3-D.
slice3D, slicecont3D: slices from a full 3-D data set.
isosurf3D: isosurfaces from a full 3-D data set as triangles.
voxel3D: isosurfaces from a full 3-D data set as points.
surf3D, spheresurf3D: 3-D shapes or surfaces.
arrows3D: arrows in 3-D.
segments3D: line segments in 3-D.
polygon3D: 3-D polygons.
box3D, border3D, rect3D: boxes and rectangles in 3-D.
text3D: labels in 3-D.
Functions defined on the image function:
image2D, for an image function to visualise 2-D or 3-D data.
ImageOcean: an image of the ocean's bathymetry.
Other plotting functions:
contour2D, for a contour function to visualise 2-D data and that have a color key.
scatter2D: colored points, lines, ... in 2-D.
text2D, arrows2D, segments2D, rect2D, polygon2D for other 2D functions that have a color key.
Colors and colorkey:
colkey: adds a color legend.
jet.col, ramp.col, gg.col, alpha.col
: suitable colors, shade and lighting.
Utility functions:
mesh: to generate rectangular (x, y) or (x, y, z) meshes.
Data sets:
Oxsat: 3-D data set with the ocean's oxygen saturation values.
Hypsometry: 2-D data set with the worlds elevation and ocean's bathymetry.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 | # run all examples
## Not run:
example(persp3D)
example(surf3D)
example(slice3D)
example(scatter3D)
example(segments3D)
example(image2D)
example(image3D)
example(contour3D)
example(colkey)
example(jet.col)
example(perspbox)
example(mesh)
example(trans3D)
example(plot.plist)
example(ImageOcean)
example(Oxsat)
## End(Not run)
|
prsp3D> # save plotting parameters
prsp3D> pm <- par("mfrow")
prsp3D> ## =======================================================================
prsp3D> ## Ribbon, persp, color keys, facets
prsp3D> ## =======================================================================
prsp3D>
prsp3D> par(mfrow = c(2, 2))
prsp3D> # simple, no scaling, use breaks to set colors
prsp3D> persp3D(z = volcano, main = "volcano", clab = c("height", "m"),
prsp3D+ breaks = seq(80,200, by = 10))
prsp3D> # keep ratios between x- and y (scale = FALSE)
prsp3D> # change ratio between x- and z (expand)
prsp3D> persp3D(z = volcano, x = 1: nrow(volcano), y = 1:ncol(volcano),
prsp3D+ expand = 0.3, main = "volcano", facets = FALSE, scale = FALSE,
prsp3D+ clab = "height, m", colkey = list(side = 1, length = 0.5))
prsp3D> # ribbon, in x--direction
prsp3D> V <- volcano[, seq(1, ncol(volcano), by = 3)] # lower resolution
prsp3D> ribbon3D(z = V, colkey = list(width = 0.5, length = 0.5,
prsp3D+ cex.axis = 0.8, side = 2), clab = "m")
prsp3D> # ribbon, in y-direction
prsp3D> Vy <- volcano[seq(1, nrow(volcano), by = 3), ]
prsp3D> ribbon3D(z = Vy, expand = 0.3, space = 0.3, along = "y",
prsp3D+ colkey = list(width = 0.5, length = 0.5, cex.axis = 0.8))
prsp3D> ## =======================================================================
prsp3D> ## Several ways to visualise 3-D data
prsp3D> ## =======================================================================
prsp3D>
prsp3D> x <- seq(-pi, pi, by = 0.2)
prsp3D> y <- seq(-pi, pi, by = 0.3)
prsp3D> grid <- mesh(x, y)
prsp3D> z <- with(grid, cos(x) * sin(y))
prsp3D> par(mfrow = c(2,2))
prsp3D> persp3D(z = z, x = x, y = y)
prsp3D> persp3D(z = z, x = x, y = y, facets = FALSE, curtain = TRUE)
prsp3D> # ribbons in two directions and larger spaces
prsp3D> ribbon3D(z = z, x = x, y = y, along = "xy", space = 0.3)
prsp3D> hist3D(z = z, x = x, y = y, border = "black")
prsp3D> ## =======================================================================
prsp3D> ## Contours and images added
prsp3D> ## =======================================================================
prsp3D>
prsp3D> par(mfrow = c(2, 2))
prsp3D> x <- seq(1, nrow(volcano), by = 3)
prsp3D> y <- seq(1, ncol(volcano), by = 3)
prsp3D> Volcano <- volcano [x, y]
prsp3D> ribbon3D(z = Volcano, contour = TRUE, zlim= c(-100, 200),
prsp3D+ image = TRUE)
prsp3D> persp3D(z = Volcano, contour = TRUE, zlim= c(-200, 200), image = FALSE)
prsp3D> persp3D(z = Volcano, x = x, y = y, scale = FALSE,
prsp3D+ contour = list(nlevels = 20, col = "red"),
prsp3D+ zlim = c(-200, 200), expand = 0.2,
prsp3D+ image = list(col = grey (seq(0, 1, length.out = 100))))
prsp3D> persp3D(z = Volcano, contour = list(side = c("zmin", "z", "350")),
prsp3D+ zlim = c(-100, 400), phi = 20, image = list(side = 350))
prsp3D> ## =======================================================================
prsp3D> ## Use of inttype
prsp3D> ## =======================================================================
prsp3D>
prsp3D> par(mfrow = c(2, 2))
prsp3D> persp3D(z = Volcano, shade = 0.5, colkey = FALSE)
prsp3D> persp3D(z = Volcano, inttype = 2, shade = 0.5, colkey = FALSE)
prsp3D> x <- y <- seq(0, 2*pi, length.out = 10)
prsp3D> z <- with (mesh(x, y), cos(x) *sin(y)) + runif(100)
prsp3D> cv <- matrix(nrow = 10, 0.5*runif(100))
prsp3D> persp3D(x, y, z, colvar = cv) # takes averages of z
prsp3D> persp3D(x, y, z, colvar = cv, inttype = 2) # takes averages of colvar
prsp3D> ## =======================================================================
prsp3D> ## Use of inttype with NAs
prsp3D> ## =======================================================================
prsp3D>
prsp3D> par(mfrow = c(2, 2))
prsp3D> VV <- V2 <- volcano[10:15, 10:15]
prsp3D> V2[3:4, 3:4] <- NA
prsp3D> V2[4, 5] <- NA
prsp3D> image2D(V2, border = "black") # shows true NA region
prsp3D> # averages of V2, including NAs, NA region larger
prsp3D> persp3D(z = VV, colvar = V2, inttype = 1, theta = 0,
prsp3D+ phi = 20, border = "black", main = "inttype = 1")
prsp3D> # extension of VV; NAs unaffected
prsp3D> persp3D(z = VV, colvar = V2, inttype = 2, theta = 0,
prsp3D+ phi = 20, border = "black", main = "inttype = 2")
prsp3D> # average of V2, ignoring NA; NA region smaller
prsp3D> persp3D(z = VV, colvar = V2, inttype = 3, theta = 0,
prsp3D+ phi = 20, border = "black", main = "inttype = 3")
prsp3D> ## =======================================================================
prsp3D> ## Use of panel.first
prsp3D> ## =======================================================================
prsp3D>
prsp3D> par(mfrow = c(1, 1))
prsp3D> # A function that is called after the axes were drawn
prsp3D> panelfirst <- function(trans) {
prsp3D+ zticks <- seq(100, 180, by = 20)
prsp3D+ len <- length(zticks)
prsp3D+ XY0 <- trans3D(x = rep(1, len), y = rep(1, len), z = zticks,
prsp3D+ pmat = trans)
prsp3D+ XY1 <- trans3D(x = rep(1, len), y = rep(61, len), z = zticks,
prsp3D+ pmat = trans)
prsp3D+ segments(XY0$x, XY0$y, XY1$x, XY1$y, lty = 2)
prsp3D+
prsp3D+ rm <- rowMeans(volcano)
prsp3D+ XY <- trans3D(x = 1:87, y = rep(ncol(volcano), 87),
prsp3D+ z = rm, pmat = trans)
prsp3D+ lines(XY, col = "blue", lwd = 2)
prsp3D+ }
prsp3D> persp3D(z = volcano, x = 1:87, y = 1: 61, scale = FALSE, theta = 10,
prsp3D+ expand = 0.2, panel.first = panelfirst, colkey = FALSE)
prsp3D> ## =======================================================================
prsp3D> ## with / without colvar / facets
prsp3D> ## =======================================================================
prsp3D>
prsp3D> par(mfrow = c(2, 2))
prsp3D> persp3D(z = volcano, shade = 0.3, col = gg.col(100))
prsp3D> # shiny colors - set lphi for more brightness
prsp3D> persp3D(z = volcano, lighting = TRUE, lphi = 90)
prsp3D> persp3D(z = volcano, col = "lightblue", colvar = NULL,
prsp3D+ shade = 0.3, bty = "b2")
prsp3D> # this also works:
prsp3D> # persp3D(z = volcano, col = "grey", shade = 0.3)
prsp3D>
prsp3D> # tilted x- and y-coordinates of 'volcano'
prsp3D> volcx <- matrix(nrow = 87, ncol = 61, data = 1:87)
prsp3D> volcx <- volcx + matrix(nrow = 87, ncol = 61,
prsp3D+ byrow = TRUE, data = seq(0., 15, length.out = 61))
prsp3D> volcy <- matrix(ncol = 87, nrow = 61, data = 1:61)
prsp3D> volcy <- t(volcy + matrix(ncol = 87, nrow = 61,
prsp3D+ byrow = TRUE, data = seq(0., 15, length.out = 87)))
prsp3D> persp3D(volcano, x = volcx, y = volcy, phi = 80)
prsp3D> ## =======================================================================
prsp3D> ## Several persps on one plot
prsp3D> ## =======================================================================
prsp3D>
prsp3D> par(mfrow = c(1, 1))
prsp3D> clim <- range(volcano)
prsp3D> persp3D(z = volcano, zlim = c(100, 600), clim = clim,
prsp3D+ box = FALSE, plot = FALSE)
prsp3D> persp3D(z = volcano + 200, clim = clim, colvar = volcano,
prsp3D+ add = TRUE, colkey = FALSE, plot = FALSE)
prsp3D> persp3D(z = volcano + 400, clim = clim, colvar = volcano,
prsp3D+ add = TRUE, colkey = FALSE) # plot = TRUE by default
prsp3D> ## =======================================================================
prsp3D> ## hist3D
prsp3D> ## =======================================================================
prsp3D>
prsp3D> par(mfrow = c(2, 2))
prsp3D> VV <- volcano[seq(1, 87, 15), seq(1, 61, 15)]
prsp3D> hist3D(z = VV, scale = FALSE, expand = 0.01, border = "black")
prsp3D> # transparent colors
prsp3D> hist3D(z = VV, scale = FALSE, expand = 0.01,
prsp3D+ alpha = 0.5, opaque.top = TRUE, border = "black")
prsp3D> hist3D(z = VV, scale = FALSE, expand = 0.01, facets = FALSE, lwd = 2)
prsp3D> hist3D(z = VV, scale = FALSE, expand = 0.01, facets = NA)
prsp3D> ## =======================================================================
prsp3D> ## hist3D and ribbon3D with greyish background, rotated, rescaled,...
prsp3D> ## =======================================================================
prsp3D>
prsp3D> par(mfrow = c(2, 2))
prsp3D> hist3D(z = VV, scale = FALSE, expand = 0.01, bty = "g", phi = 20,
prsp3D+ col = "#0072B2", border = "black", shade = 0.2, ltheta = 90,
prsp3D+ space = 0.3, ticktype = "detailed", d = 2)
prsp3D> # extending the ranges
prsp3D> plotdev(xlim = c(-0.2, 1.2), ylim = c(-0.2, 1.2), theta = 45)
prsp3D> ribbon3D(z = VV, scale = FALSE, expand = 0.01, bty = "g", phi = 20,
prsp3D+ col = "lightblue", border = "black", shade = 0.2, ltheta = 90,
prsp3D+ space = 0.3, ticktype = "detailed", d = 2, curtain = TRUE)
prsp3D> ribbon3D(z = VV, scale = FALSE, expand = 0.01, bty = "g", phi = 20, zlim = c(95,183),
prsp3D+ col = "lightblue", lighting = TRUE, ltheta = 50, along = "y",
prsp3D+ space = 0.7, ticktype = "detailed", d = 2, curtain = TRUE)
prsp3D> ## =======================================================================
prsp3D> ## hist3D for a 1-D data set
prsp3D> ## =======================================================================
prsp3D>
prsp3D> par(mfrow = c(2, 1))
prsp3D> x <- rchisq(1000, df = 4)
prsp3D> hs <- hist(x, breaks = 15)
prsp3D> hist3D(x = hs$mids, y = 1, z = matrix(ncol = 1, data = hs$density),
prsp3D+ bty = "g", ylim = c(0., 2.0), scale = FALSE, expand = 20,
prsp3D+ border = "black", col = "white", shade = 0.3, space = 0.1,
prsp3D+ theta = 20, phi = 20, main = "3-D perspective")
prsp3D> # reset plotting parameters
prsp3D> par(mfrow = pm)
surf3D> # save plotting parameters
surf3D> pm <- par("mfrow")
surf3D> pmar <- par("mar")
surf3D> par(mar = c(1, 1, 1, 1))
surf3D> ## =======================================================================
surf3D> ## A three-dimensional shape
surf3D> ## (ala http://docs.enthought.com/mayavi/mayavi/mlab.html)
surf3D> ## =======================================================================
surf3D>
surf3D> par(mfrow = c(2, 2))
surf3D> # create grid matrices
surf3D> X <- seq(0, pi, length.out = 50)
surf3D> Y <- seq(0, 2*pi, length.out = 50)
surf3D> M <- mesh(X, Y)
surf3D> phi <- M$x
surf3D> theta <- M$y
surf3D> # x, y and z grids
surf3D> r <- sin(4*phi)^3 + cos(2*phi)^3 + sin(6*theta)^2 + cos(6*theta)^4
surf3D> x <- r * sin(phi) * cos(theta)
surf3D> y <- r * cos(phi)
surf3D> z <- r * sin(phi) * sin(theta)
surf3D> # full colored image
surf3D> surf3D(x, y, z, colvar = y, colkey = FALSE, shade = 0.5,
surf3D+ box = FALSE, theta = 60)
surf3D> # same, but just facets
surf3D> surf3D(x, y, z, colvar = y, colkey = FALSE, box = FALSE,
surf3D+ theta = 60, facets = FALSE)
surf3D> # with colors and border, AND increasing the size
surf3D> # (by reducing the x- y and z- ranges
surf3D> surf3D(x, y, z, colvar = y, colkey = FALSE, box = FALSE,
surf3D+ theta = 60, border = "black", xlim = range(x)*0.8,
surf3D+ ylim = range(y)*0.8, zlim = range(z)*0.8)
surf3D> # Now with one color and shading
surf3D> surf3D(x, y, z, box = FALSE,
surf3D+ theta = 60, col = "lightblue", shade = 0.9)
surf3D> ## Not run:
surf3D> ##D # rotation
surf3D> ##D for (angle in seq(0, 360, by = 10))
surf3D> ##D plotdev(theta = angle)
surf3D> ##D
surf3D> ## End(Not run)
surf3D>
surf3D> ## =======================================================================
surf3D> ## Several other shapes
surf3D> ## http://xahlee.info/surface/gallery.html
surf3D> ## =======================================================================
surf3D>
surf3D> par(mfrow = c(2, 2))
surf3D> # Shape 1
surf3D> M <- mesh(seq(0, 6*pi, length.out = 50),
surf3D+ seq(pi/3, pi, length.out = 50))
surf3D> u <- M$x ; v <- M$y
surf3D> x <- u/2 * sin(v) * cos(u)
surf3D> y <- u/2 * sin(v) * sin(u)
surf3D> z <- u/2 * cos(v)
surf3D> surf3D(x, y, z, colvar = z, colkey = FALSE, box = FALSE, phi = 50)
surf3D> # Shape 2: add border
surf3D> M <- mesh(seq(0, 2*pi, length.out = 50),
surf3D+ seq(0, 2*pi, length.out = 50))
surf3D> u <- M$x ; v <- M$y
surf3D> x <- sin(u)
surf3D> y <- sin(v)
surf3D> z <- sin(u + v)
surf3D> surf3D(x, y, z, colvar = z, border = "black",
surf3D+ colkey = FALSE)
surf3D> # shape 3: uses same mesh, other perspective (d >1)
surf3D> x <- (3 + cos(v/2)*sin(u) - sin(v/2)*sin(2*u))*cos(v)
surf3D> y <- (3 + cos(v/2)*sin(u) - sin(v/2)*sin(2*u))*sin(v)
surf3D> z <- sin(v/2)*sin(u) + cos(v/2)*sin(2*u)
surf3D> surf3D(x, y, z, colvar = z, colkey = FALSE, d = 2, facets = FALSE)
surf3D> # shape 4: more complex colvar
surf3D> M <- mesh(seq(-13.2, 13.2, length.out = 50),
surf3D+ seq(-37.4, 37.4, length.out = 50))
surf3D> u <- M$x ; v <- M$y
surf3D> b <- 0.4; r <- 1 - b^2; w <- sqrt(r)
surf3D> D <- b*((w*cosh(b*u))^2 + (b*sin(w*v))^2)
surf3D> x <- -u + (2*r*cosh(b*u)*sinh(b*u)) / D
surf3D> y <- (2*w*cosh(b*u)*(-(w*cos(v)*cos(w*v)) - sin(v)*sin(w*v))) / D
surf3D> z <- (2*w*cosh(b*u)*(-(w*sin(v)*cos(w*v)) + cos(v)*sin(w*v))) / D
surf3D> surf3D(x, y, z, colvar = sqrt(x + 8.3), colkey = FALSE,
surf3D+ theta = 10, border = "black", box = FALSE)
surf3D> box()
surf3D> ## =======================================================================
surf3D> ## A sphere, with box type with grid lines
surf3D> ## =======================================================================
surf3D>
surf3D> par(mar = c(2, 2, 2, 2))
surf3D> par(mfrow = c(1, 1))
surf3D> M <- mesh(seq(0, 2*pi, length.out = 50),
surf3D+ seq(0, pi, length.out = 50))
surf3D> u <- M$x ; v <- M$y
surf3D> x <- cos(u)*sin(v)
surf3D> y <- sin(u)*sin(v)
surf3D> z <- cos(v)
surf3D> colvar <- sin(u*6) * sin(v*6)
surf3D> surf3D(y, x, z, colvar = colvar, phi = 0, bty = "b2",
surf3D+ lighting = TRUE, ltheta = 40)
surf3D> ## =======================================================================
surf3D> ## Function spheresurf3D
surf3D> ## =======================================================================
surf3D>
surf3D> par(mfrow = c(2, 2))
surf3D> spheresurf3D()
surf3D> # true ranges are [-1, 1]; set limits to [-0.8, 0.8] to make larger plots
surf3D> lim <- c(-0.8, 0.8)
surf3D> spheresurf3D(colkey = FALSE, xlim = lim, ylim = lim, zlim = lim)
surf3D> spheresurf3D(bty = "b", ticktype = "detailed", phi = 50)
surf3D> spheresurf3D(colvar = matrix(nrow = 30, data = runif(900)))
surf3D> ## =======================================================================
surf3D> ## Images on a sphere
surf3D> ## =======================================================================
surf3D>
surf3D> par(mfrow = c(1, 1), mar = c(1, 1, 1, 3))
surf3D> AA <- Hypsometry$z; AA[AA<=0] <- NA
surf3D> lim <- c(-0.8, 0.8)
surf3D> # log transformation of color variable
surf3D> spheresurf3D(AA, NAcol = "black", theta = 90, phi = 30, box = FALSE,
surf3D+ xlim = lim, ylim = lim, zlim = lim, log = "c")
surf3D> # restore plotting parameters
surf3D> par(mfrow = pm)
surf3D> par(mar = pmar)
slic3D> # save plotting parameters
slic3D> pm <- par("mfrow")
slic3D> pmar <- par("mar")
slic3D> ## =======================================================================
slic3D> ## Simple slice3D examples
slic3D> ## =======================================================================
slic3D>
slic3D> par(mfrow = c(2, 2))
slic3D> x <- y <- z <- seq(-1, 1, by = 0.1)
slic3D> grid <- mesh(x, y, z)
slic3D> colvar <- with(grid, x*exp(-x^2 - y^2 - z^2))
slic3D> # default is just the panels
slic3D> slice3D (x, y, z, colvar = colvar, theta = 60)
slic3D> # contour slices
slic3D> slicecont3D (x, y, z, ys = seq(-1, 1, by = 0.5), colvar = colvar,
slic3D+ theta = 60, border = "black")
slic3D> slice3D (x, y, z, xs = c(-1, -0.5, 0.5), ys = c(-1, 0, 1),
slic3D+ zs = c(-1, 0), colvar = colvar,
slic3D+ theta = 60, phi = 40)
slic3D> ## =======================================================================
slic3D> ## coloring on a surface
slic3D> ## =======================================================================
slic3D>
slic3D> XY <- mesh(x, y)
slic3D> ZZ <- XY$x*XY$y
slic3D> slice3D (x, y, z, xs = XY$x, ys = XY$y, zs = ZZ, colvar = colvar,
slic3D+ lighting = TRUE, lphi = 90, ltheta = 0)
slic3D> ## =======================================================================
slic3D> ## Specifying transparent colors
slic3D> ## =======================================================================
slic3D>
slic3D> par(mfrow = c(1, 1))
slic3D> x <- y <- z <- seq(-4, 4, by = 0.2)
slic3D> M <- mesh(x, y, z)
slic3D> R <- with (M, sqrt(x^2 + y^2 + z^2))
slic3D> p <- sin(2*R) /(R+1e-3)
slic3D> ## Not run:
slic3D> ##D # This is very slow - alpha = 0.5 makes it transparent
slic3D> ##D
slic3D> ##D slice3D(x, y, z, colvar = p, col = jet.col(alpha = 0.5),
slic3D> ##D xs = 0, ys = c(-4, 0, 4), zs = NULL, d = 2)
slic3D> ## End(Not run)
slic3D>
slic3D> slice3D(x, y, z, colvar = p, d = 2, theta = 60, border = "black",
slic3D+ xs = c(-4, 0), ys = c(-4, 0, 4), zs = c(-4, 0))
slic3D> ## =======================================================================
slic3D> ## A section along a transect
slic3D> ## =======================================================================
slic3D>
slic3D> data(Oxsat)
slic3D> Ox <- Oxsat$val[, Oxsat$lat > - 5 & Oxsat$lat < 5, ]
slic3D> slice3D(x = Oxsat$lon, z = -Oxsat$depth, y = 1:5, colvar = Ox,
slic3D+ ys = 1:5, zs = NULL, NAcol = "black",
slic3D+ expand = 0.4, theta = 45, phi = 45)
slic3D> ## =======================================================================
slic3D> ## isosurf3D example - rather slow
slic3D> ## =======================================================================
slic3D>
slic3D> par(mfrow = c(2, 2), mar = c(2, 2, 2, 2))
slic3D> x <- y <- z <- seq(-2, 2, length.out = 15)
slic3D> xyz <- mesh(x, y, z)
slic3D> F <- with(xyz, log(x^2 + y^2 + z^2 +
slic3D+ 10*(x^2 + y^2) * (y^2 + z^2) ^2))
slic3D> # use shading for level = 1 - show triangulation with border
slic3D> isosurf3D(x, y, z, F, level = 1, shade = 0.9,
slic3D+ col = "yellow", border = "orange")
slic3D> # lighting for level - 2
slic3D> isosurf3D(x, y, z, F, level = 2, lighting = TRUE,
slic3D+ lphi = 0, ltheta = 0, col = "blue", shade = NA)
slic3D> # three levels, transparency added
slic3D> isosurf3D(x, y, z, F, level = seq(0, 4, by = 2),
slic3D+ col = c("red", "blue", "yellow"),
slic3D+ clab = "F", alpha = 0.2, theta = 0, lighting = TRUE)
slic3D> # transparency can also be added afterwards with plotdev()
slic3D> ## Not run:
slic3D> ##D isosurf3D(x, y, z, F, level = seq(0, 4, by = 2),
slic3D> ##D col = c("red", "blue", "yellow"),
slic3D> ##D shade = NA, plot = FALSE, clab = "F")
slic3D> ##D plotdev(lighting = TRUE, alpha = 0.2, theta = 0)
slic3D> ## End(Not run)
slic3D> # use of creatisosurf
slic3D> iso <- createisosurf(x, y, z, F, level = 2)
slic3D> head(iso)
x y z
[1,] 0.0000000 -0.1179802 -2.000000
[2,] -0.1204879 0.0000000 -2.000000
[3,] 0.0000000 -0.2073773 -1.714286
[4,] 0.0000000 -0.2073773 -1.714286
[5,] -0.1204879 0.0000000 -2.000000
[6,] -0.2138894 0.0000000 -1.714286
slic3D> triangle3D(iso, col = "green", shade = 0.3)
slic3D> ## Not run:
slic3D> ##D # higher resolution
slic3D> ##D x <- y <- z <- seq(-2, 2, length.out = 50)
slic3D> ##D xyz <- mesh(x, y, z)
slic3D> ##D F <- with(xyz, log(x^2 + y^2 + z^2 +
slic3D> ##D 10*(x^2 + y^2) * (y^2 + z^2) ^2))
slic3D> ##D
slic3D> ##D # three levels
slic3D> ##D isosurf3D(x, y, z, F, level = seq(0, 4, by = 2),
slic3D> ##D col = c("red", "blue", "yellow"),
slic3D> ##D shade = NA, plot = FALSE, clab = "F")
slic3D> ##D plotdev(lighting = TRUE, alpha = 0.2, theta = 0)
slic3D> ## End(Not run)
slic3D>
slic3D> ## =======================================================================
slic3D> ## voxel3D example
slic3D> ## =======================================================================
slic3D>
slic3D> par(mfrow = c(2, 2), mar = c(2, 2, 2, 2))
slic3D> # fast but needs high resolution grid
slic3D> x <- y <- z <- seq(-2, 2, length.out = 70)
slic3D> xyz <- mesh(x, y, z)
slic3D> F <- with(xyz, log(x^2 + y^2 + z^2 +
slic3D+ 10*(x^2 + y^2) * (y^2 + z^2) ^2))
slic3D> voxel3D(x, y, z, F, level = 4, pch = ".", cex = 5)
slic3D> ## =======================================================================
slic3D> ## rotation
slic3D> ## =======================================================================
slic3D>
slic3D> plotdev(theta = 45, phi = 0)
slic3D> plotdev(theta = 90, phi = 10)
slic3D> # same using createvoxel - more flexible for coloring
slic3D> vox <- createvoxel(x, y, z, F, level = 4)
slic3D> scatter3D(vox$x, vox$y, vox$z, colvar = vox$y,
slic3D+ bty = "g", colkey = FALSE)
slic3D> ## =======================================================================
slic3D> ## voxel3D to show hypox sites
slic3D> ## =======================================================================
slic3D>
slic3D> par(mfrow = c(1, 1), mar = c(2, 2, 2, 2))
slic3D> Hypox <- createvoxel(Oxsat$lon, Oxsat$lat, Oxsat$depth[1:19],
slic3D+ Oxsat$val[,,1:19], level = 40, operator = "<")
slic3D> panel <- function(pmat) { # an image at the bottom
slic3D+ Nx <- length(Oxsat$lon)
slic3D+ Ny <- length(Oxsat$lat)
slic3D+ M <- mesh(Oxsat$lon, Oxsat$lat)
slic3D+ xy <- trans3D(pmat = pmat, x = as.vector(M$x), y = as.vector(M$y),
slic3D+ z = rep(-1000, length.out = Nx*Ny))
slic3D+ x <- matrix(nrow = Nx, ncol = Ny, data = xy$x)
slic3D+ y <- matrix(nrow = Nx, ncol = Ny, data = xy$y)
slic3D+ Bat <- Oxsat$val[,,1]; Bat[!is.na(Bat)] <- 1
slic3D+ image2D(x = x, y = y, z = Bat, NAcol = "black", col = "grey",
slic3D+ add = TRUE, colkey = FALSE)
slic3D+ }
slic3D> scatter3D(Hypox$x, Hypox$y, -Hypox$z, colvar = Hypox$cv,
slic3D+ panel.first = panel, pch = ".", bty = "b",
slic3D+ theta = 30, phi = 20, ticktype = "detailed",
slic3D+ zlim = c(-1000,0), xlim = range(Oxsat$lon),
slic3D+ ylim = range(Oxsat$lat) )
slic3D> # restore plotting parameters
slic3D> par(mfrow = pm)
slic3D> par(mar = pmar)
sctt3D> # save plotting parameters
sctt3D> pm <- par("mfrow")
sctt3D> ## =======================================================================
sctt3D> ## A sphere
sctt3D> ## =======================================================================
sctt3D>
sctt3D> par(mfrow = c(1, 1))
sctt3D> M <- mesh(seq(0, 2*pi, length.out = 100),
sctt3D+ seq(0, pi, length.out = 100))
sctt3D> u <- M$x ; v <- M$y
sctt3D> x <- cos(u)*sin(v)
sctt3D> y <- sin(u)*sin(v)
sctt3D> z <- cos(v)
sctt3D> # full panels of box are drawn (bty = "f")
sctt3D> scatter3D(x, y, z, pch = ".", col = "red",
sctt3D+ bty = "f", cex = 2, colkey = FALSE)
sctt3D> ## =======================================================================
sctt3D> ## Different types
sctt3D> ## =======================================================================
sctt3D>
sctt3D> par (mfrow = c(2, 2))
sctt3D> z <- seq(0, 10, 0.2)
sctt3D> x <- cos(z)
sctt3D> y <- sin(z)*z
sctt3D> # greyish background for the boxtype (bty = "g")
sctt3D> scatter3D(x, y, z, phi = 0, bty = "g",
sctt3D+ pch = 20, cex = 2, ticktype = "detailed")
sctt3D> # add another point
sctt3D> scatter3D(x = 0, y = 0, z = 0, add = TRUE, colkey = FALSE,
sctt3D+ pch = 18, cex = 3, col = "black")
sctt3D> # add text
sctt3D> text3D(x = cos(1:10), y = (sin(1:10)*(1:10) - 1),
sctt3D+ z = 1:10, colkey = FALSE, add = TRUE,
sctt3D+ labels = LETTERS[1:10], col = c("black", "red"))
sctt3D> # line plot
sctt3D> scatter3D(x, y, z, phi = 0, bty = "g", type = "l",
sctt3D+ ticktype = "detailed", lwd = 4)
sctt3D> # points and lines
sctt3D> scatter3D(x, y, z, phi = 0, bty = "g", type = "b",
sctt3D+ ticktype = "detailed", pch = 20,
sctt3D+ cex = c(0.5, 1, 1.5))
sctt3D> # vertical lines
sctt3D> scatter3D(x, y, z, phi = 0, bty = "g", type = "h",
sctt3D+ ticktype = "detailed")
sctt3D> ## =======================================================================
sctt3D> ## With confidence interval
sctt3D> ## =======================================================================
sctt3D>
sctt3D> x <- runif(20)
sctt3D> y <- runif(20)
sctt3D> z <- runif(20)
sctt3D> par(mfrow = c(1, 1))
sctt3D> CI <- list(z = matrix(nrow = length(x),
sctt3D+ data = rep(0.05, 2*length(x))))
sctt3D> # greyish background for the boxtype (bty = "g")
sctt3D> scatter3D(x, y, z, phi = 0, bty = "g", CI = CI,
sctt3D+ col = gg.col(100), pch = 18, cex = 2, ticktype = "detailed",
sctt3D+ xlim = c(0, 1), ylim = c(0, 1), zlim = c(0, 1))
sctt3D> # add new set of points
sctt3D> x <- runif(20)
sctt3D> y <- runif(20)
sctt3D> z <- runif(20)
sctt3D> CI2 <- list(x = matrix(nrow = length(x),
sctt3D+ data = rep(0.05, 2*length(x))),
sctt3D+ z = matrix(nrow = length(x),
sctt3D+ data = rep(0.05, 2*length(x))))
sctt3D> scatter3D(x, y, z, CI = CI2, add = TRUE, col = "red", pch = 16)
sctt3D> ## =======================================================================
sctt3D> ## With a surface
sctt3D> ## =======================================================================
sctt3D>
sctt3D> par(mfrow = c(1, 1))
sctt3D> # surface = volcano
sctt3D> M <- mesh(1:nrow(volcano), 1:ncol(volcano))
sctt3D> # 100 points above volcano
sctt3D> N <- 100
sctt3D> xs <- runif(N) * 87
sctt3D> ys <- runif(N) * 61
sctt3D> zs <- runif(N)*50 + 154
sctt3D> # scatter + surface
sctt3D> scatter3D(xs, ys, zs, ticktype = "detailed", pch = 16,
sctt3D+ bty = "f", xlim = c(1, 87), ylim = c(1,61), zlim = c(94, 215),
sctt3D+ surf = list(x = M$x, y = M$y, z = volcano,
sctt3D+ NAcol = "grey", shade = 0.1))
sctt3D> ## =======================================================================
sctt3D> ## A surface and CI
sctt3D> ## =======================================================================
sctt3D>
sctt3D> par(mfrow = c(1, 1))
sctt3D> M <- mesh(seq(0, 2*pi, length = 30), (1:30)/100)
sctt3D> z <- with (M, sin(x) + y)
sctt3D> # points 'sampled'
sctt3D> N <- 30
sctt3D> xs <- runif(N) * 2*pi
sctt3D> ys <- runif(N) * 0.3
sctt3D> zs <- sin(xs) + ys + rnorm(N)*0.3
sctt3D> CI <- list(z = matrix(nrow = length(xs),
sctt3D+ data = rep(0.3, 2*length(xs))),
sctt3D+ lwd = 3)
sctt3D> # facets = NA makes a transparent surface; borders are black
sctt3D> scatter3D(xs, ys, zs, ticktype = "detailed", pch = 16,
sctt3D+ xlim = c(0, 2*pi), ylim = c(0, 0.3), zlim = c(-1.5, 1.5),
sctt3D+ CI = CI, theta = 20, phi = 30, cex = 2,
sctt3D+ surf = list(x = M$x, y = M$y, z = z, border = "black", facets = NA)
sctt3D+ )
sctt3D> ## =======================================================================
sctt3D> ## droplines till the fitted surface
sctt3D> ## =======================================================================
sctt3D>
sctt3D> with (mtcars, {
sctt3D+
sctt3D+ # linear regression
sctt3D+ fit <- lm(mpg ~ wt + disp)
sctt3D+
sctt3D+ # predict values on regular xy grid
sctt3D+ wt.pred <- seq(1.5, 5.5, length.out = 30)
sctt3D+ disp.pred <- seq(71, 472, length.out = 30)
sctt3D+ xy <- expand.grid(wt = wt.pred,
sctt3D+ disp = disp.pred)
sctt3D+
sctt3D+ mpg.pred <- matrix (nrow = 30, ncol = 30,
sctt3D+ data = predict(fit, newdata = data.frame(xy),
sctt3D+ interval = "prediction"))
sctt3D+
sctt3D+ # fitted points for droplines to surface
sctt3D+ fitpoints <- predict(fit)
sctt3D+
sctt3D+ scatter3D(z = mpg, x = wt, y = disp, pch = 18, cex = 2,
sctt3D+ theta = 20, phi = 20, ticktype = "detailed",
sctt3D+ xlab = "wt", ylab = "disp", zlab = "mpg",
sctt3D+ surf = list(x = wt.pred, y = disp.pred, z = mpg.pred,
sctt3D+ facets = NA, fit = fitpoints),
sctt3D+ main = "mtcars")
sctt3D+
sctt3D+ })
sctt3D> ## =======================================================================
sctt3D> ## Two ways to make a scatter 3D of quakes data set
sctt3D> ## =======================================================================
sctt3D>
sctt3D> par(mfrow = c(1, 1))
sctt3D> # first way, use vertical spikes (type = "h")
sctt3D> with(quakes, scatter3D(x = long, y = lat, z = -depth, colvar = mag,
sctt3D+ pch = 16, cex = 1.5, xlab = "longitude", ylab = "latitude",
sctt3D+ zlab = "depth, km", clab = c("Richter","Magnitude"),
sctt3D+ main = "Earthquakes off Fiji", ticktype = "detailed",
sctt3D+ type = "h", theta = 10, d = 2,
sctt3D+ colkey = list(length = 0.5, width = 0.5, cex.clab = 0.75))
sctt3D+ )
sctt3D> # second way: add dots on bottom and left panel
sctt3D> # before the scatters are drawn,
sctt3D> # add small dots on basal plane and on the depth plane
sctt3D> panelfirst <- function(pmat) {
sctt3D+ zmin <- min(-quakes$depth)
sctt3D+ XY <- trans3D(quakes$long, quakes$lat,
sctt3D+ z = rep(zmin, nrow(quakes)), pmat = pmat)
sctt3D+ scatter2D(XY$x, XY$y, colvar = quakes$mag, pch = ".",
sctt3D+ cex = 2, add = TRUE, colkey = FALSE)
sctt3D+
sctt3D+ xmin <- min(quakes$long)
sctt3D+ XY <- trans3D(x = rep(xmin, nrow(quakes)), y = quakes$lat,
sctt3D+ z = -quakes$depth, pmat = pmat)
sctt3D+ scatter2D(XY$x, XY$y, colvar = quakes$mag, pch = ".",
sctt3D+ cex = 2, add = TRUE, colkey = FALSE)
sctt3D+ }
sctt3D> with(quakes, scatter3D(x = long, y = lat, z = -depth, colvar = mag,
sctt3D+ pch = 16, cex = 1.5, xlab = "longitude", ylab = "latitude",
sctt3D+ zlab = "depth, km", clab = c("Richter","Magnitude"),
sctt3D+ main = "Earthquakes off Fiji", ticktype = "detailed",
sctt3D+ panel.first = panelfirst, theta = 10, d = 2,
sctt3D+ colkey = list(length = 0.5, width = 0.5, cex.clab = 0.75))
sctt3D+ )
sctt3D> ## =======================================================================
sctt3D> ## text3D and scatter3D
sctt3D> ## =======================================================================
sctt3D>
sctt3D> with(USArrests, text3D(Murder, Assault, Rape,
sctt3D+ colvar = UrbanPop, col = gg.col(100), theta = 60, phi = 20,
sctt3D+ xlab = "Murder", ylab = "Assault", zlab = "Rape",
sctt3D+ main = "USA arrests",
sctt3D+ labels = rownames(USArrests), cex = 0.6,
sctt3D+ bty = "g", ticktype = "detailed", d = 2,
sctt3D+ clab = c("Urban","Pop"), adj = 0.5, font = 2))
sctt3D> with(USArrests, scatter3D(Murder, Assault, Rape - 1,
sctt3D+ colvar = UrbanPop, col = gg.col(100),
sctt3D+ type = "h", pch = ".", add = TRUE))
sctt3D> ## =======================================================================
sctt3D> ## zoom near origin
sctt3D> ## =======================================================================
sctt3D>
sctt3D> # display axis ranges
sctt3D> getplist()[c("xlim","ylim","zlim")]
$xlim
[1] 0.8 17.4
$ylim
[1] 45 337
$zlim
[1] 7.3 46.0
sctt3D> # choose suitable ranges
sctt3D> plotdev(xlim = c(0, 10), ylim = c(40, 150),
sctt3D+ zlim = c(7, 25))
sctt3D> ## =======================================================================
sctt3D> ## text3D to label x- and y axis
sctt3D> ## =======================================================================
sctt3D>
sctt3D> par(mfrow = c(1, 1))
sctt3D> hist3D (x = 1:5, y = 1:4, z = VADeaths,
sctt3D+ bty = "g", phi = 20, theta = -60,
sctt3D+ xlab = "", ylab = "", zlab = "", main = "VADeaths",
sctt3D+ col = "#0072B2", border = "black", shade = 0.8,
sctt3D+ ticktype = "detailed", space = 0.15, d = 2, cex.axis = 1e-9)
sctt3D> text3D(x = 1:5, y = rep(0.5, 5), z = rep(3, 5),
sctt3D+ labels = rownames(VADeaths),
sctt3D+ add = TRUE, adj = 0)
sctt3D> text3D(x = rep(1, 4), y = 1:4, z = rep(0, 4),
sctt3D+ labels = colnames(VADeaths),
sctt3D+ add = TRUE, adj = 1)
sctt3D> ## =======================================================================
sctt3D> ## Scatter2D; bty can also be set = to one of the perspbox alernatives
sctt3D> ## =======================================================================
sctt3D>
sctt3D> par(mfrow = c(2, 2))
sctt3D> x <- seq(0, 2*pi, length.out = 30)
sctt3D> scatter2D(x, sin(x), colvar = cos(x), pch = 16,
sctt3D+ ylab = "sin", clab = "cos", cex = 1.5)
sctt3D> # other box types:
sctt3D> scatter2D(x, sin(x), colvar = cos(x), type = "l", lwd = 4, bty = "g")
sctt3D> scatter2D(x, sin(x), colvar = cos(x), type = "b", lwd = 2, bty = "b2")
sctt3D> # transparent colors and spikes
sctt3D> scatter2D(x, sin(x), colvar = cos(x), type = "h", lwd = 4, alpha = 0.5)
sctt3D> ## =======================================================================
sctt3D> ## mesh examples and scatter2D
sctt3D> ## =======================================================================
sctt3D>
sctt3D> par(mfrow = c(1, 2))
sctt3D> x <- seq(-1, 1, by = 0.1)
sctt3D> y <- seq(-2, 2, by = 0.2)
sctt3D> grid <- mesh(x, y)
sctt3D> z <- with(grid, cos(x) * sin(y))
sctt3D> image2D(z, x = x, y = y)
sctt3D> points(grid)
sctt3D> scatter2D(grid$x, grid$y, colvar = z, pch = 20, cex = 2)
sctt3D> ## =======================================================================
sctt3D> ## scatter plot with confidence intervals
sctt3D> ## =======================================================================
sctt3D>
sctt3D> par(mfrow = c(2, 2))
sctt3D> x <- sort(rnorm(10))
sctt3D> y <- runif(10)
sctt3D> cv <- sqrt(x^2 + y^2)
sctt3D> CI <- list(lwd = 2)
sctt3D> CI$x <- matrix (nrow = length(x), data = c(rep(0.25, 2*length(x))))
sctt3D> scatter2D(x, y, colvar = cv, pch = 16, cex = 2, CI = CI)
sctt3D> scatter2D(x, y, colvar = cv, pch = 16, cex = 2, CI = CI, type = "b")
sctt3D> CI$y <- matrix (nrow = length(x), data = c(rep(0.05, 2*length(x))))
sctt3D> CI$col <- "black"
sctt3D> scatter2D(x, y, colvar = cv, pch = 16, cex = 2, CI = CI)
sctt3D> CI$y[c(2,4,8,10), ] <- NA # Some points have no CI
sctt3D> CI$x[c(2,4,8,10), ] <- NA # Some points have no CI
sctt3D> CI$alen <- 0.02 # increase arrow head
sctt3D> scatter2D(x, y, colvar = cv, pch = 16, cex = 2, CI = CI)
sctt3D> ## =======================================================================
sctt3D> ## Scatter on an image
sctt3D> ## =======================================================================
sctt3D>
sctt3D> par(mfrow = c(1, 1))
sctt3D> # image of oxygen saturation
sctt3D> oxlim <- range(Oxsat$val[,,1], na.rm = TRUE)
sctt3D> image2D(z = Oxsat$val[,,1], x = Oxsat$lon, y = Oxsat$lat,
sctt3D+ contour = TRUE,
sctt3D+ xlab = "longitude", ylab = "latitude",
sctt3D+ main = "Oxygen saturation", clim = oxlim, clab = "%")
sctt3D> # (imaginary) measurements at 5 sites
sctt3D> lon <- c( 11.2, 6.0, 0.9, -4, -8.8)
sctt3D> lat <- c(-19.7,-14.45,-9.1,-3.8, -1.5)
sctt3D> O2sat <- c( 90, 95, 92, 85, 100)
sctt3D> # add to image; use same zrange; avoid adding a color key
sctt3D> scatter2D(colvar = O2sat, x = lon, y = lat, clim = oxlim, pch = 16,
sctt3D+ add = TRUE, cex = 2, colkey = FALSE)
sctt3D> ## =======================================================================
sctt3D> ## Scatter on a contourplot
sctt3D> ## =======================================================================
sctt3D>
sctt3D> par(mfrow = c(1, 1))
sctt3D> # room for colorkey by setting colkey = list(plot = FALSE)
sctt3D>
sctt3D> # contour plot of the ocean's bathymetry
sctt3D> Depth <- Hypsometry$z
sctt3D> Depth[Depth > 0] <- NA
sctt3D> contour2D(z = Depth, x = Hypsometry$x, y = Hypsometry$y,
sctt3D+ xlab = "longitude", ylab = "latitude",
sctt3D+ col = "black", NAcol = "grey", levels = seq(-6000, 0, by = 2000),
sctt3D+ main = "Oxygen saturation along ship track",
sctt3D+ colkey = list(plot = FALSE))
sctt3D> # add data to image; with a color key
sctt3D> scatter2D(colvar = O2sat, x = lon, y = lat, pch = 16,
sctt3D+ add = TRUE, cex = 2, clab = "%")
sctt3D> ## =======================================================================
sctt3D> ## scatter2D for time-series plots
sctt3D> ## =======================================================================
sctt3D>
sctt3D> # Plotting sunspot 'anomalies'
sctt3D> sunspot <- data.frame(year = time(sunspot.month),
sctt3D+ anom = sunspot.month - mean(sunspot.month))
sctt3D> # long-term moving average of anomaly
sctt3D> ff <- 100
sctt3D> sunspot$ma <- filter(sunspot$anom, rep(1/ff, ff), sides = 2)
sctt3D> with (sunspot, lines2D(year, anom,
sctt3D+ colvar = anom > 0,
sctt3D+ col = c("pink", "lightblue"),
sctt3D+ main = "sunspot anomaly", type = "h",
sctt3D+ colkey = FALSE, las = 1, xlab = "year", ylab = ""))
sctt3D> lines2D(sunspot$year, sunspot$ma, add = TRUE)
sctt3D> # The same
sctt3D> #with (sunspot, plot(year, anom,
sctt3D> # col = c("pink", "lightblue")[(anom > 0) + 1],
sctt3D> # main = "sunspot", type = "h", las = 1))
sctt3D>
sctt3D> # but this does not work due to NAs...
sctt3D> # lines(sunspot$year, sunspot$ma)
sctt3D>
sctt3D> ## =======================================================================
sctt3D> ## text2D
sctt3D> ## =======================================================================
sctt3D>
sctt3D> with(USArrests, text2D(x = Murder, y = Assault + 5, colvar = Rape,
sctt3D+ xlab = "Murder", ylab = "Assault", clab = "Rape",
sctt3D+ main = "USA arrests", labels = rownames(USArrests), cex = 0.6,
sctt3D+ adj = 0.5, font = 2))
sctt3D> with(USArrests, scatter2D(x = Murder, y = Assault, colvar = Rape,
sctt3D+ pch = 16, add = TRUE, colkey = FALSE))
sctt3D> # reset plotting parameters
sctt3D> par(mfrow = pm)
Warning message:
In cbind(Col[-1], Col[-len]) :
number of rows of result is not a multiple of vector length (arg 1)
sgmn3D> # save plotting parameters
sgmn3D> pm <- par("mfrow")
sgmn3D> ## ========================================================================
sgmn3D> ## arrows, points, segments, box
sgmn3D> ## ========================================================================
sgmn3D>
sgmn3D> # Create a grid of x, y, and z values
sgmn3D> xx <- yy <- seq(-0.8, 0.8, by = 0.2)
sgmn3D> zz <- seq(-0.8, 0.8, by = 0.8)
sgmn3D> M <- mesh(xx, yy, zz)
sgmn3D> x0 <- M$x; y0 <- M$y; z0 <- M$z
sgmn3D> x1 <- x0 + 0.1
sgmn3D> Col <- c("red", "blue", "green")
sgmn3D> arrows3D(x0, y0, z0, x1 = x1, colvar = z0, lwd = 2,
sgmn3D+ d = 2, clab = "z-value", col = Col, length = 0.1,
sgmn3D+ xlim = c(-0.8, 0.8), ylim = c(-0.8, 0.8),
sgmn3D+ main = "arrows3D, points3D, segments3D, border3D")
sgmn3D> # add starting point of arrows
sgmn3D> points3D(x0, y0, z0, add = TRUE, colvar = z0,
sgmn3D+ colkey = FALSE, pch = ".", cex = 3, col = Col)
sgmn3D> # use segments to add section
sgmn3D> x0 <- c(-0.8, 0.8, 0.8, -0.8)
sgmn3D> x1 <- c( 0.8, 0.8, -0.8, -0.8)
sgmn3D> y0 <- c(-0.8, -0.8, 0.8, -0.8)
sgmn3D> y1 <- c(-0.8, 0.8, 0.8, 0.8)
sgmn3D> z0 <- c(0., 0., 0., 0.)
sgmn3D> segments3D(x0, y0, z0, x1, y1, z1 = z0,
sgmn3D+ add = TRUE, col = "black", lwd = 2)
sgmn3D> # add a box
sgmn3D> border3D(-0.8, -0.8, -0.8, 0.8, 0.8, 0.8,
sgmn3D+ col = "orange", add = TRUE, lwd = 3)
sgmn3D> ## ========================================================================
sgmn3D> ## boxes, cubes
sgmn3D> ## ========================================================================
sgmn3D>
sgmn3D> # borders are boxes without facets
sgmn3D> border3D(x0 = seq(-0.8, -0.1, by = 0.1),
sgmn3D+ y0 = seq(-0.8, -0.1, by = 0.1),
sgmn3D+ z0 = seq(-0.8, -0.1, by = 0.1),
sgmn3D+ x1 = seq(0.8, 0.1, by = -0.1),
sgmn3D+ y1 = seq(0.8, 0.1, by = -0.1),
sgmn3D+ z1 = seq(0.8, 0.1, by = -0.1),
sgmn3D+ col = gg.col(8), lty = 2,
sgmn3D+ lwd = c(1, 4), phi = 20, main = "border3D")
sgmn3D> box3D(x0 = -0.8, y0 = -0.8, z0 = -0.8,
sgmn3D+ x1 = 0.8, y1 = 0.8, z1 = 0.8,
sgmn3D+ border = "black", lwd = 2,
sgmn3D+ col = gg.col(1, alpha = 0.8),
sgmn3D+ main = "box3D")
sgmn3D> box3D(x0 = seq(-0.8, -0.1, by = 0.1),
sgmn3D+ y0 = seq(-0.8, -0.1, by = 0.1),
sgmn3D+ z0 = seq(-0.8, -0.1, by = 0.1),
sgmn3D+ x1 = seq(0.8, 0.1, by = -0.1),
sgmn3D+ y1 = seq(0.8, 0.1, by = -0.1),
sgmn3D+ z1 = seq(0.8, 0.1, by = -0.1),
sgmn3D+ col = rainbow(n = 8, alpha = 0.1),
sgmn3D+ border = "black", lwd = 2, phi = 20)
sgmn3D> # here the perspective does not always work
sgmn3D> # use alpha.col to set the transparency of a vector of colors
sgmn3D> box3D(x0 = runif(3), y0 = runif(3), z0 = runif(3),
sgmn3D+ x1 = runif(3), y1 = runif(3), z1 = runif(3),
sgmn3D+ col = c("red", "lightblue", "orange"), alpha = 0.5,
sgmn3D+ border = "black", lwd = 2)
sgmn3D> ## ========================================================================
sgmn3D> ## rectangles
sgmn3D> ## ========================================================================
sgmn3D> # at constant 'z'
sgmn3D> rect3D(x0 = seq(-0.8, -0.1, by = 0.1),
sgmn3D+ y0 = seq(-0.8, -0.1, by = 0.1),
sgmn3D+ z0 = seq(-0.8, -0.1, by = 0.1),
sgmn3D+ x1 = seq(0.8, 0.1, by = -0.1),
sgmn3D+ y1 = seq(0.8, 0.1, by = -0.1),
sgmn3D+ col = gg.col(8), border = "black",
sgmn3D+ bty = "g", lwd = 2, phi = 20, main = "rect3D")
sgmn3D> # constant y and with transparent facets
sgmn3D> rect3D(x0 = 0, y0 = 0, z0 = 0, x1 = 1, z1 = 5,
sgmn3D+ ylim = c(0, 1), facets = NA, border = "red",
sgmn3D+ bty = "g", lwd = 2, phi = 20)
sgmn3D> # add rect at constant z, with colored facet
sgmn3D> rect3D(x0 = 0, y0 = 0, z0 = 0, x1 = 1, y1 = 1,
sgmn3D+ border = "red", add = TRUE)
sgmn3D> ## ========================================================================
sgmn3D> ## arrows added to a persp plot
sgmn3D> ## ========================================================================
sgmn3D>
sgmn3D> x <- y <- seq(-10, 10, length = 30)
sgmn3D> z <- outer(x, y, FUN = function(x,y) x^2 + y^2)
sgmn3D> persp3D(x, y, z, theta = 30, phi = 30,
sgmn3D+ col = "lightblue", ltheta = 120, shade = 0.75,
sgmn3D+ ticktype = "detailed", xlab = "X",
sgmn3D+ ylab = "Y", zlab = "x^2+y^2" )
sgmn3D> # Points where to put the arrows
sgmn3D> x <- y <- seq(-10, 10, len = 6)
sgmn3D> X0 <- outer(x, y, FUN = function (x,y) x)
sgmn3D> Y0 <- outer(x, y, FUN = function (x,y) y)
sgmn3D> Z0 <- outer(x, y, FUN = function (x,y) x^2 + y^2)
sgmn3D> X1 <- X0 + 1
sgmn3D> Y1 <- Y0 + 1
sgmn3D> Z1 <- Z0 + 10
sgmn3D> arrows3D(X0, Y0, Z0, X1, Y1, Z1, lwd = 2,
sgmn3D+ add = TRUE, type = "curved", col = "red")
sgmn3D> segments3D(X0, Y0, Z0, X0, Y0, rep(0, length(X0)), lwd = 2,
sgmn3D+ add = TRUE, col = "green")
sgmn3D> ## ========================================================================
sgmn3D> ## polygon3D
sgmn3D> ## ========================================================================
sgmn3D>
sgmn3D> x <- runif(10)
sgmn3D> y <- runif(10)
sgmn3D> z <- runif(10)
sgmn3D> polygon3D(x, y, z)
sgmn3D> # several polygons, separated by NAs
sgmn3D> x <- runif(39)
sgmn3D> y <- runif(39)
sgmn3D> z <- runif(39)
sgmn3D> ii <- seq(4, 36, by = 4)
sgmn3D> x[ii] <- y[ii] <- z[ii] <- NA
sgmn3D> # transparent colors (alpha)
sgmn3D> polygon3D(x, y, z, border = "black", lwd = 3,
sgmn3D+ col = gg.col(length(ii) + 1, alpha = 0.8),
sgmn3D+ main = "polygon3D")
sgmn3D> ## ========================================================================
sgmn3D> ## 2D examples, with color key
sgmn3D> ## ========================================================================
sgmn3D>
sgmn3D> arrows2D(x0 = runif(10), y0 = runif(10),
sgmn3D+ x1 = runif(10), y1 = runif(10), colvar = 1:10,
sgmn3D+ code = 3, main = "arrows2D, segments2D")
sgmn3D> segments2D(x0 = runif(10), y0 = runif(10),
sgmn3D+ x1 = runif(10), y1 = runif(10), colvar = 1:10,
sgmn3D+ lwd = 2, add = TRUE, colkey = FALSE)
sgmn3D> # transparency
sgmn3D> rect2D(x0 = runif(10), y0 = runif(10),
sgmn3D+ x1 = runif(10), y1 = runif(10), colvar = 1:10,
sgmn3D+ alpha = 0.4, lwd = 2, main = "rect2D")
sgmn3D> ## ========================================================================
sgmn3D> ## polygon2D
sgmn3D> ## ========================================================================
sgmn3D>
sgmn3D> x <- runif(10)
sgmn3D> y <- runif(10)
sgmn3D> polygon2D(x, y) # same as polygon
sgmn3D> # several polygons, separated by NAs
sgmn3D> x <- runif(59)
sgmn3D> y <- runif(59)
sgmn3D> ii <- seq(5, 55, by = 5)
sgmn3D> x[ii] <- y[ii] <- NA
sgmn3D> # transparent colors (alpha)
sgmn3D> polygon2D(x, y, border = "black", lwd = 3,
sgmn3D+ colvar = 1:(length(ii) + 1),
sgmn3D+ col = gg.col(), alpha = 0.2,
sgmn3D+ main = "polygon2D")
sgmn3D> # restore plotting parameters
sgmn3D> par(mfrow = pm)
imag2D> # save plotting parameters
imag2D> pm <- par("mfrow")
imag2D> ## =======================================================================
imag2D> ## Difference between x or y a vector/matrix and rasterImage
imag2D> ## =======================================================================
imag2D>
imag2D> par(mfrow = c(2, 2))
imag2D> x <- y <- 1:3
imag2D> z <- matrix (nrow = 3, ncol = 3, data = 1:9)
imag2D> image2D(z, x, y, border = "black")
imag2D> image2D(z, x, y, rasterImage = TRUE, border = "black")
imag2D> image2D(z, x = matrix(nrow = 3, ncol = 3, data = x),
imag2D+ y, border = "black")
imag2D> image2D(z, x, y, border = "black", theta = 45)
imag2D> ## =======================================================================
imag2D> ## shading, light, adding contours, points and lines
imag2D> ## =======================================================================
imag2D>
imag2D> par(mfrow = c(2, 2))
imag2D> nr <- nrow(volcano)
imag2D> nc <- ncol(volcano)
imag2D> image2D(volcano, x = 1:nr, y = 1:nc, lighting = TRUE,
imag2D+ main = "volcano", clab = "height, m")
imag2D> abline(v = seq(10, 80, by = 10))
imag2D> abline(h = seq(10, 60, by = 10))
imag2D> points(50, 30, pch = 3, cex = 5, lwd = 3, col = "white")
imag2D> image2D(z = volcano, x = 1:nr, y = 1:nc, lwd = 2, shade = 0.2,
imag2D+ main = "volcano", clab = "height, m")
imag2D> image2D(volcano, x = 1:nr, y = 1:nc, contour = TRUE, shade = 0.5, lphi = 0,
imag2D+ col = "lightblue", main = "volcano")
imag2D> breaks <- seq(90, 200, by = 10)
imag2D> image2D(volcano, x = 1:nr, y = 1:nc, col = jet.col(length(breaks)-1),
imag2D+ main = "volcano", clab = "height, m", breaks = breaks)
imag2D> ## =======================================================================
imag2D> ## Contour plots
imag2D> ## =======================================================================
imag2D>
imag2D> par(mfrow = c(2, 2))
imag2D> V <- volcano - 150
imag2D> # default, no color key
imag2D> contour2D(z = V, colkey = FALSE, lwd = 2)
imag2D> # imposed levels
imag2D> contour2D(z = V, lwd = 2, levels = seq(-40, 40, by = 20))
imag2D> # negative levels dashed
imag2D> contour2D(z = V, col = "black", lwd = 2,
imag2D+ levels = seq(0, 40, by = 20))
imag2D> contour2D(z = V, col = "black", lwd = 2, lty = 2,
imag2D+ levels = seq(-40, -20, by = 20), add = TRUE)
imag2D> # no labels, imposed number of levels, colorkey
imag2D> contour2D(z = V, lwd = 2, nlevels = 20, drawlabels = FALSE,
imag2D+ colkey = list(at = seq(-40, 40, by = 20)))
imag2D> ## =======================================================================
imag2D> ## A large data set, input is an array
imag2D> ## =======================================================================
imag2D>
imag2D> par(mfrow = c(1, 1))
imag2D> image2D(z = Oxsat$val[, , 1], x = Oxsat$lon, y = Oxsat$lat,
imag2D+ main = "surface oxygen saturation data 2005", NAcol = "black",
imag2D+ clab = c("","","%"))
imag2D> # images at first 9 depths - use subset to select them
imag2D> image2D(z = Oxsat$val, subset = 1:9,
imag2D+ x = Oxsat$lon, y = Oxsat$lat,
imag2D+ margin = c(1, 2), NAcol = "black",
imag2D+ xlab = "longitude", ylab = "latitude",
imag2D+ zlim = c(0, 115),
imag2D+ main = paste("depth ", Oxsat$depth[1:9], " m"),
imag2D+ mfrow = c(3, 3))
imag2D> # images at latitude - depth section - increase resolution
imag2D> z <- Oxsat$val[, Oxsat$lat > - 5 & Oxsat$lat < 5, ]
imag2D> image2D(z = z, x = Oxsat$lon, y = Oxsat$depth,
imag2D+ margin = c(1, 3), NAcol = "black",
imag2D+ resfac = 3, ylim = c(5000, 0))
imag2D> # show position of transects
imag2D> image2D(z = Oxsat$val[ , ,1],
imag2D+ x = Oxsat$lon, y = Oxsat$lat,
imag2D+ NAcol = "black")
imag2D> abline(h = Oxsat$lat[Oxsat$lat > - 5 & Oxsat$lat < 5])
imag2D> ## =======================================================================
imag2D> ## Image of a list of matrices
imag2D> ## =======================================================================
imag2D>
imag2D> listvolcano <- list(volcano = volcano, logvolcano = log(volcano))
imag2D> image2D(listvolcano, x = 1:nr, y = 1:nc, contour = TRUE,
imag2D+ main = c("volcano", "log(volcano)"),
imag2D+ clab = list("height, m", "log(m)"),
imag2D+ zlim = list(c(80, 200), c(4.4, 5.5)))
imag2D> ## =======================================================================
imag2D> ## Image of a list of arrays
imag2D> ## =======================================================================
imag2D>
imag2D> ## Not run:
imag2D> ##D # crude conversion from oxsat to oxygen
imag2D> ##D listoxygen <- list(Oxsat$val, Oxsat$val/100 * 360)
imag2D> ##D
imag2D> ##D image2D(z = listoxygen,
imag2D> ##D x = Oxsat$lon, y = Oxsat$lat,
imag2D> ##D margin = c(1, 2), NAcol = "black",
imag2D> ##D main = c("Oxygen saturation ", " Oxygen concentration"),
imag2D> ##D mtext = paste("depth ", Oxsat$depth, " m")
imag2D> ##D )
imag2D> ## End(Not run)
imag2D>
imag2D> ## =======================================================================
imag2D> ## 'x', 'y' and 'z' are matrices
imag2D> ## =======================================================================
imag2D>
imag2D> par(mfrow = c(2, 1))
imag2D> # tilted x- and y-coordinates of 'volcano'
imag2D> volcx <- matrix(nrow = 87, ncol = 61, data = 1:87)
imag2D> volcx <- volcx + matrix(nrow = 87, ncol = 61,
imag2D+ byrow = TRUE, data = seq(0., 15, length.out = 61))
imag2D> volcy <- matrix(ncol = 87, nrow = 61, data = 1:61)
imag2D> volcy <- t(volcy + matrix(ncol = 87, nrow = 61,
imag2D+ byrow = TRUE, data = seq(0., 25, length.out = 87)))
imag2D> image2D(volcano, x = volcx, y = volcy)
imag2D> # x and y can also be of dimension dim(z)+1:
imag2D> ## Not run:
imag2D> ##D # tilted x- and y-coordinates of 'volcano'
imag2D> ##D volcx <- matrix(nrow = 88, ncol = 62, data = 1:88)
imag2D> ##D volcx <- volcx + matrix(nrow = 88, ncol = 62,
imag2D> ##D byrow = TRUE, data = seq(0., 15, length.out = 62))
imag2D> ##D
imag2D> ##D volcy <- matrix(ncol = 88, nrow = 62, data = 1:62)
imag2D> ##D volcy <- t(volcy + matrix(ncol = 88, nrow = 62,
imag2D> ##D byrow = TRUE, data = seq(0., 25, length.out = 88)))
imag2D> ##D
imag2D> ##D image2D(volcano, x = volcx, y = volcy)
imag2D> ## End(Not run)
imag2D>
imag2D> # use of panel function
imag2D> image2D(volcano, x = volcx, y = volcy, NAcol = "black",
imag2D+ panel.first = substitute(box(col = "lightgrey", lwd = 30)))
imag2D> ## =======================================================================
imag2D> ## Image with NAs and logs
imag2D> ## =======================================================================
imag2D>
imag2D> par(mfrow = c(2, 2))
imag2D> # normal volcano
imag2D> image2D(volcano, clab = c("height", "m"))
imag2D> # logarithmic z-axis
imag2D> image2D(volcano, log = "z", clab = c("height", "m"),
imag2D+ main = "log='z'")
imag2D> # Including NAs
imag2D> VOLC <- volcano - 110
imag2D> VOLC [VOLC <= 0] <- NA
imag2D> image2D(VOLC, main = "including NAs and rescaled")
imag2D> # both
imag2D> image2D(VOLC, NAcol = "black", log = "z", zlim = c(1, 100),
imag2D+ main = "NAs and log = 'z'")
imag2D> ## =======================================================================
imag2D> ## Image with contour specification (alpha sets the transparency)
imag2D> ## =======================================================================
imag2D>
imag2D> par(mfrow = c(1, 1))
imag2D> image2D(volcano, shade = 0.2, rasterImage = TRUE,
imag2D+ contour = list(col = "white", labcex = 0.8, lwd = 3, alpha = 0.5))
imag2D> # same:
imag2D> ## Not run:
imag2D> ##D image2D(z = volcano, shade = 0.2, rasterImage = TRUE)
imag2D> ##D contour2D(z = volcano, col = "white", labcex = 0.8,
imag2D> ##D lwd = 3, alpha = 0.5, add = TRUE)
imag2D> ## End(Not run)
imag2D> # reset plotting parameters
imag2D> par(mfrow = pm)
imag3D> # save plotting parameters
imag3D> pm <- par("mfrow")
imag3D> ## =======================================================================
imag3D> ## images in x, y, z plane
imag3D> ## =======================================================================
imag3D>
imag3D> par(mfrow = c(2, 2))
imag3D> # images in x, y, z plane
imag3D> # We use colkey = list(plot = FALSE) to create room for a color key
imag3D> image3D(y = seq(0, 1, 0.1), z = seq(0, 1, 0.1), x = 0.5,
imag3D+ col = "blue", xlim = c(0,1), colkey = list(plot = FALSE))
imag3D> image3D(x = seq(0, 1, 0.1), z = seq(0, 1, 0.1), y = 0.5,
imag3D+ add = TRUE, col = "red", alpha = 0.2) # alpha makes it transparent
imag3D> image3D(x = seq(0, 1, 0.1), y = seq(0, 1, 0.1), z = 0.5,
imag3D+ add = TRUE, col = "green")
imag3D> colkey(col = c("green", "red", "blue"), clim = c(0.5, 3.5),
imag3D+ at = 1:3, labels = c("z", "y", "x"), add = TRUE)
imag3D> #
imag3D> image3D(z = 100, colvar = volcano, zlim = c(0, 150),
imag3D+ clab = c("height", "m"))
imag3D> #
imag3D> image3D( x = 0.5, colvar = volcano, xlim = c(0, 1),
imag3D+ ylim = c(0, 1), zlim = c(0, 1))
imag3D> image3D( y = 0.5, colvar = volcano, add = TRUE)
imag3D> #
imag3D> image3D( z = 1, colvar = volcano,
imag3D+ x = seq(0, 1, length.out = nrow(volcano)),
imag3D+ y = seq(0, 1, length.out = ncol(volcano)),
imag3D+ xlim = c(0, 2), ylim = c(0, 2), zlim = c(0, 2))
imag3D> image3D(y = 2, colvar = volcano, add = TRUE,
imag3D+ shade = 0.2,
imag3D+ x = seq(0, 1, length.out = nrow(volcano)),
imag3D+ z = seq(1, 2, length.out = ncol(volcano)))
imag3D> image3D(x = 2, colvar = NULL, col = "orange", add = TRUE,
imag3D+ y = seq(0, 1, length.out = nrow(volcano)),
imag3D+ z = seq(1, 2, length.out = ncol(volcano)))
imag3D> # reset plotting parameters
imag3D> par(mfrow = pm)
cntr3D> # save plotting parameters
cntr3D> pm <- par("mfrow")
cntr3D> ## =======================================================================
cntr3D> ## Contours
cntr3D> ## =======================================================================
cntr3D> par (mfrow = c(2, 2))
cntr3D> r <- 1:nrow(volcano)
cntr3D> c <- 1:ncol(volcano)
cntr3D> contour3D(x = r, y = c, z = 100, colvar = volcano, zlim = c(0, 150),
cntr3D+ clab = c("height", "m"))
cntr3D> contour3D(x = 100, y = r, z = c, colvar = volcano, clab = c("height", "m"))
cntr3D> contour3D(z = volcano, colvar = volcano, lwd = 2,
cntr3D+ nlevels = 20, clab = c("height", "m"), colkey = FALSE)
cntr3D> contour3D(y = volcano, colvar = volcano, lwd = 2,
cntr3D+ nlevels = 10, clab = c("height", "m"))
cntr3D> ## =======================================================================
cntr3D> ## Composite images and contours in 3D
cntr3D> ## =======================================================================
cntr3D> persp3D(z = volcano, zlim = c(90, 300), col = "white",
cntr3D+ shade = 0.1, d = 2, plot = FALSE)
cntr3D> contour3D(z = volcano, colvar = volcano, lwd = 2, add = TRUE,
cntr3D+ nlevels = 20, clab = c("height", "m"), plot = FALSE,
cntr3D+ colkey = list(at = seq(90, 190, length.out = 5)))
cntr3D> contour3D(z = 300, colvar = volcano, lwd = 2, col = "grey",
cntr3D+ add = TRUE, nlevels = 5)
cntr3D> ## =======================================================================
cntr3D> ## the viewing depth of contours (dDepth)
cntr3D> ## =======================================================================
cntr3D>
cntr3D> # too low
cntr3D> persp3D(z = volcano, col = "white", shade = 0.1, plot = FALSE)
cntr3D> contour3D(z = volcano, colvar = volcano, lwd = 2,
cntr3D+ add = TRUE, dDepth = 0, col = "black")
cntr3D> # default
cntr3D> persp3D(z = volcano, col = "white", shade = 0.1, plot = FALSE)
cntr3D> contour3D(z = volcano, colvar = volcano, lwd = 2,
cntr3D+ add = TRUE, dDepth = 0.1, col = "black")
cntr3D> # too high
cntr3D> persp3D(z = volcano, col = "white", shade = 0.1, plot = FALSE)
cntr3D> contour3D(z = volcano, colvar = volcano, lwd = 1,
cntr3D+ add = TRUE, dDepth = 0.5, col = "black")
cntr3D> # reset plotting parameters
cntr3D> par(mfrow = pm)
colkey> # save plotting parameters
colkey> pm <- par(mfrow = c(2, 2))
colkey> pmar <- par(mar = c(5.1, 4.1, 4.1, 2.1))
colkey> ## =======================================================================
colkey> ## colorkey as argument of a plot3D function
colkey> ## =======================================================================
colkey> # default, colkey = NULL: adds colkey because multiple colors
colkey> image2D(z = volcano)
colkey> # default, colkey = NULL: no colkey because only one color
colkey> image2D(z = volcano, col = "grey", shade = 0.2, contour = TRUE)
colkey> # colkey = FALSE: no color key, no extra space foreseen
colkey> image2D(z = volcano, colkey = FALSE)
colkey> # colkey = list(plot = FALSE): no color key, extra space foreseen
colkey> image2D(z = volcano, colkey = list(plot = FALSE, side = 3))
colkey> colkey (side = 3, add = TRUE, clim = range(volcano))
colkey> ## =======================================================================
colkey> ## colorkey in new plot
colkey> ## =======================================================================
colkey>
colkey> colkey(side = 1, clim = c(0, 1), add = FALSE, clab = "z",
colkey+ col.clab = "red", adj.clab = 0)
colkey> colkey(side = 2, clim = c(0, 1), clab = "z", length = 0.5, width = 0.5)
colkey> colkey(side = 3, clim = c(0, 1), lwd = 3, clab = c("a","b","c","d"),
colkey+ line.clab = 5)
colkey> colkey(side = 4, clim = c(1e-6, 1), clog = TRUE,
colkey+ clab = "a very long title in bold and close to the key",
colkey+ line.clab = 1, side.clab = 2, font.clab = 2)
colkey> ## =======================================================================
colkey> ## colorkey added to existing plot
colkey> ## =======================================================================
colkey>
colkey> par(mfrow = c(1, 1))
colkey> image2D(volcano, xlab = "", clab = "m",
colkey+ colkey = list(side = 1, length = 0.5, width = 0.5,
colkey+ line.clab = 1))
colkey> colkey(side = 3, clim = range(volcano), add = TRUE)
colkey> # 'dist' to put colkey within the image
colkey> # 'shift' to position colkey to the right or upward
colkey> par(mfrow = c(1, 1))
colkey> image2D(volcano, colkey = FALSE)
colkey> colkey(clim = range(volcano), dist = -0.15, shift = 0.2,
colkey+ side = 3, add = TRUE, clab = "key 1", col.clab = "white",
colkey+ length = 0.5, width = 0.5, col.axis = "white",
colkey+ col.ticks = "white", cex.axis = 0.8)
colkey> colkey(clim = range(volcano), dist = -0.1, shift = -0.2,
colkey+ side = 4, add = TRUE, clab = "key 2", col.clab = "white",
colkey+ length = 0.3, width = 0.5, col.axis = "white",
colkey+ col.ticks = "white", col.box = "red", cex.axis = 0.8)
colkey> colkey(clim = range(volcano), dist = -0.3,
colkey+ side = 1, add = TRUE, clab = "key 3", col.clab = "white",
colkey+ length = 0.3, width = 0.5, col.axis = "white",
colkey+ col.ticks = "white", at = c(100, 140, 180),
colkey+ labels = c("a", "b", "c"), font = 2)
colkey> colkey(clim = range(volcano), dist = -0.3, shift = -0.2,
colkey+ side = 2, add = TRUE, clab = "key 4", col.clab = "white",
colkey+ length = 0.3, width = 0.5, col.axis = "white",
colkey+ col.ticks = "white", col.box = "red", cex.axis = 0.8,
colkey+ las = 3)
colkey> ## =======================================================================
colkey> ## colorkey in other plots
colkey> ## =======================================================================
colkey>
colkey> par(mfrow = c(1, 1))
colkey> par(mar = par("mar") + c(0, 0, -2, 0))
colkey> image2D(volcano, clab = "height, m",
colkey+ colkey = list(dist = -0.15, shift = 0.2,
colkey+ side = 3, length = 0.5, width = 0.5, line.clab = 2.5,
colkey+ cex.clab = 2, col.clab = "white", col.axis = "white",
colkey+ col.ticks = "white", cex.axis = 0.8))
colkey> ## =======================================================================
colkey> ## Several color keys in composite plot
colkey> ## =======================================================================
colkey>
colkey> persp3D(z = volcano, zlim = c(-60, 200), phi = 20, bty = "b",
colkey+ colkey = list(length = 0.2, width = 0.4, shift = 0.15,
colkey+ cex.axis = 0.8, cex.clab = 0.85), lighting = TRUE, lphi = 90,
colkey+ clab = c("height","m"), plot = FALSE)
colkey> # create gradient in x-direction
colkey> Vx <- volcano[-1, ] - volcano[-nrow(volcano), ]
colkey> # add as image with own color key, at bottom
colkey> image3D(z = -60, colvar = Vx/10, add = TRUE,
colkey+ colkey = list(length = 0.2, width = 0.4, shift = -0.15,
colkey+ cex.axis = 0.8, cex.clab = 0.85),
colkey+ clab = c("gradient","m/m"), plot = TRUE)
colkey> ## =======================================================================
colkey> ## categorical colors; use addlines = TRUE to separate colors
colkey> ## =======================================================================
colkey>
colkey> with(iris, scatter3D(x = Sepal.Length, y = Sepal.Width,
colkey+ z = Petal.Length, colvar = as.integer(Species),
colkey+ col = c("orange", "green", "lightblue"), pch = 16, cex = 2,
colkey+ clim = c(1, 3), ticktype = "detailed", phi = 20,
colkey+ xlab = "Sepal Length", ylab = "Sepal Width",
colkey+ zlab = "Petal Length", main = "iris",
colkey+ colkey = list(at = c(1.33, 2, 2.66), side = 1,
colkey+ addlines = TRUE, length = 0.5, width = 0.5,
colkey+ labels = c("setosa", "versicolor", "virginica") )))
colkey> # reset plotting parameters
colkey> par(mfrow = pm)
colkey> par(mar = pmar)
jet.cl> # save plotting parameters
jet.cl> pm <- par("mfrow")
jet.cl> pmar <- par("mar")
jet.cl> ## =======================================================================
jet.cl> ## Transparency and various color schemes
jet.cl> ## =======================================================================
jet.cl>
jet.cl> par(mfrow = c(3, 3))
jet.cl> for (alph in c(0.25, 0.75))
jet.cl+ image2D(volcano, alpha = alph,
jet.cl+ main = paste("jet.col, alpha = ", alph))
jet.cl> image2D(volcano, main = "jet.col")
jet.cl> image2D(volcano, col = jet2.col(100), main = "jet2.col")
jet.cl> image2D(volcano, col = gg.col(100), main = "gg.col")
jet.cl> image2D(volcano, col = gg2.col(100), main = "gg2.col")
jet.cl> image2D(volcano, col = rainbow(100), main = "rainbow")
jet.cl> image2D(volcano, col = terrain.colors(100), main = "terrain.colors")
jet.cl> image2D(volcano, col = ramp.col(c("blue", "yellow", "green", "red")),
jet.cl+ main = "ramp.col")
jet.cl> ## =======================================================================
jet.cl> ## Shading, lighting - one color
jet.cl> ## =======================================================================
jet.cl>
jet.cl> # create grid matrices
jet.cl> X <- seq(0, pi, length.out = 50)
jet.cl> Y <- seq(0, 2*pi, length.out = 50)
jet.cl> M <- mesh(X, Y)
jet.cl> phi <- M$x
jet.cl> theta <- M$y
jet.cl> # x, y and z grids
jet.cl> x <- sin(phi) * cos(theta)
jet.cl> y <- cos(phi)
jet.cl> z <- sin(phi) * sin(theta)
jet.cl> # these are the defaults
jet.cl> p <- list(ambient = 0.3, diffuse = 0.6, specular = 1.,
jet.cl+ exponent = 20, sr = 0, alpha = 1)
jet.cl> par(mfrow = c(3, 3), mar = c(0, 0, 0, 0))
jet.cl> Col <- "red"
jet.cl> surf3D(x, y, z, box = FALSE, col = Col, shade = 0.9)
jet.cl> surf3D(x, y, z, box = FALSE, col = Col, lighting = TRUE)
jet.cl> surf3D(x, y, z, box = FALSE, col = Col, lighting = list(ambient = 0))
jet.cl> surf3D(x, y, z, box = FALSE, col = Col, lighting = list(diffuse = 0))
jet.cl> surf3D(x, y, z, box = FALSE, col = Col, lighting = list(diffuse = 1))
jet.cl> surf3D(x, y, z, box = FALSE, col = Col, lighting = list(specular = 0))
jet.cl> surf3D(x, y, z, box = FALSE, col = Col, lighting = list(exponent = 5))
jet.cl> surf3D(x, y, z, box = FALSE, col = Col, lighting = list(exponent = 50))
jet.cl> surf3D(x, y, z, box = FALSE, col = Col, lighting = list(sr = 1))
jet.cl> ## =======================================================================
jet.cl> ## Shading, lighting with default colors
jet.cl> ## =======================================================================
jet.cl>
jet.cl> x <- seq(-pi, pi, len = 100)
jet.cl> y <- seq(-pi, pi, len = 100)
jet.cl> grid <- mesh(x, y)
jet.cl> z <- with(grid, cos(x) * sin(y))
jet.cl> cv <- with(grid, -cos(y) * sin(x))
jet.cl> # lphi = 180, ltheta = -130 - good for shade
jet.cl> # lphi = 90, ltheta = 0 - good for lighting
jet.cl>
jet.cl> par(mfrow = c(2, 2))
jet.cl> persp3D(z = z, x = x, y = y, colvar = cv, zlim = c(-3, 3), colkey = FALSE)
jet.cl> persp3D(z = z, x = x, y = y, colvar = cv, zlim = c(-3, 3),
jet.cl+ lighting = TRUE, colkey = FALSE)
jet.cl> persp3D(z = z, x = x, y = y, colvar = cv, zlim = c(-3, 3),
jet.cl+ shade = 0.25, colkey = FALSE)
jet.cl> persp3D(z = z, x = x, y = y, colvar = cv, zlim = c(-3, 3),
jet.cl+ lighting = TRUE, lphi = 90, ltheta = 0, colkey = FALSE)
jet.cl> ## =======================================================================
jet.cl> ## transparency of a vector of colors
jet.cl> ## =======================================================================
jet.cl>
jet.cl> par(mfrow = c(1, 1))
jet.cl> x <- runif(19)
jet.cl> y <- runif(19)
jet.cl> z <- runif(19)
jet.cl> # split into 5 sections (polygons)
jet.cl> ii <- seq(4, 19, by = 4)
jet.cl> x[ii] <- y[ii] <- z[ii] <- NA
jet.cl> polygon3D(x, y, z, border = "black", lwd = 2,
jet.cl+ col = alpha.col(c("red", "lightblue", "yellow", "green", "black"),
jet.cl+ alpha = 0.4))
jet.cl> # the same, now passing alpha as an argument to polygon3D:
jet.cl> ## Not run:
jet.cl> ##D polygon3D(x, y, z, border = "black", lwd = 2,
jet.cl> ##D col = c("red", "lightblue", "yellow", "green", "black"),
jet.cl> ##D alpha = 0.4)
jet.cl> ## End(Not run)
jet.cl> # reset plotting parameters
jet.cl> par(mfrow = pm)
jet.cl> par(mar = pmar)
prspbx> # save plotting parameters
prspbx> pm <- par("mfrow")
prspbx> pmar <- par("mar")
prspbx> ## ========================================================================
prspbx> ## The 4 predefined box types
prspbx> ## ========================================================================
prspbx>
prspbx> par(mfrow = c(2, 2), mar = c(1, 1, 1, 1))
prspbx> # box type with only backward panels
prspbx> perspbox(z = volcano, bty = "b", ticktype = "detailed", d = 2,
prspbx+ main = "bty = 'b'")
prspbx> # box as in 'persp'
prspbx> perspbox(z = volcano, bty = "f", ticktype = "detailed",
prspbx+ d = 2, main = "bty = 'f'")
prspbx> # back panels with gridlines, detailed axes
prspbx> perspbox(z = volcano, bty = "b2", ticktype = "detailed",
prspbx+ d = 2, main = "bty = 'b2'")
prspbx> # ggplot-type, simple axes
prspbx> perspbox(z = volcano, bty = "g",
prspbx+ d = 2, main = "bty = 'g'")
prspbx> ## ========================================================================
prspbx> ## A user-defined box
prspbx> ## ========================================================================
prspbx>
prspbx> par(mfrow = c(1, 1))
prspbx> perspbox(z = diag(2), bty = "u", ticktype = "detailed",
prspbx+ col.panel = "gold", col.axis = "white",
prspbx+ scale = FALSE, expand = 0.4,
prspbx+ col.grid = "grey", main = "user-defined")
prspbx> # restore plotting parameters
prspbx> par(mfrow = pm)
prspbx> par(mar = pmar)
mesh> ## ========================================================================
mesh> ## 2-D mesh
mesh> ## ========================================================================
mesh>
mesh> x <- c(-1 , 0, 1)
mesh> y <- 1 : 4
mesh> # 2-D mesh
mesh> (M <- mesh(x, y))
$x
[,1] [,2] [,3] [,4]
[1,] -1 -1 -1 -1
[2,] 0 0 0 0
[3,] 1 1 1 1
$y
[,1] [,2] [,3] [,4]
[1,] 1 2 3 4
[2,] 1 2 3 4
[3,] 1 2 3 4
mesh> # calculate with this mesh
mesh> V <- with (M, x/2 * sin(y))
mesh> # same as:
mesh> V2 <- outer(x, y, FUN = function(x, y) x/2*sin(y))
mesh> ## ========================================================================
mesh> ## 3-D mesh
mesh> ## ========================================================================
mesh>
mesh> x <- y <- z <- c(-1 , 0, 1)
mesh> # 3-D mesh
mesh> (M <- mesh(x, y, z))
$x
, , 1
[,1] [,2] [,3]
[1,] -1 -1 -1
[2,] 0 0 0
[3,] 1 1 1
, , 2
[,1] [,2] [,3]
[1,] -1 -1 -1
[2,] 0 0 0
[3,] 1 1 1
, , 3
[,1] [,2] [,3]
[1,] -1 -1 -1
[2,] 0 0 0
[3,] 1 1 1
$y
, , 1
[,1] [,2] [,3]
[1,] -1 0 1
[2,] -1 0 1
[3,] -1 0 1
, , 2
[,1] [,2] [,3]
[1,] -1 0 1
[2,] -1 0 1
[3,] -1 0 1
, , 3
[,1] [,2] [,3]
[1,] -1 0 1
[2,] -1 0 1
[3,] -1 0 1
$z
, , 1
[,1] [,2] [,3]
[1,] -1 -1 -1
[2,] -1 -1 -1
[3,] -1 -1 -1
, , 2
[,1] [,2] [,3]
[1,] 0 0 0
[2,] 0 0 0
[3,] 0 0 0
, , 3
[,1] [,2] [,3]
[1,] 1 1 1
[2,] 1 1 1
[3,] 1 1 1
mesh> # calculate with 3-D mesh
mesh> V <- with (M, x/2 * sin(y) *sqrt(z+2))
mesh> # plot result
mesh> scatter3D(M$x, M$y, M$z, V, pch = "+", cex = 3, colkey = FALSE)
trns3D> ## ========================================================================
trns3D> ## 3-D mesh
trns3D> ## ========================================================================
trns3D>
trns3D> x <- y <- z <- c(-1 , 0, 1)
trns3D> # plot a 3-D mesh
trns3D> (M <- mesh(x, y, z))
$x
, , 1
[,1] [,2] [,3]
[1,] -1 -1 -1
[2,] 0 0 0
[3,] 1 1 1
, , 2
[,1] [,2] [,3]
[1,] -1 -1 -1
[2,] 0 0 0
[3,] 1 1 1
, , 3
[,1] [,2] [,3]
[1,] -1 -1 -1
[2,] 0 0 0
[3,] 1 1 1
$y
, , 1
[,1] [,2] [,3]
[1,] -1 0 1
[2,] -1 0 1
[3,] -1 0 1
, , 2
[,1] [,2] [,3]
[1,] -1 0 1
[2,] -1 0 1
[3,] -1 0 1
, , 3
[,1] [,2] [,3]
[1,] -1 0 1
[2,] -1 0 1
[3,] -1 0 1
$z
, , 1
[,1] [,2] [,3]
[1,] -1 -1 -1
[2,] -1 -1 -1
[3,] -1 -1 -1
, , 2
[,1] [,2] [,3]
[1,] 0 0 0
[2,] 0 0 0
[3,] 0 0 0
, , 3
[,1] [,2] [,3]
[1,] 1 1 1
[2,] 1 1 1
[3,] 1 1 1
trns3D> # plot result
trns3D> pmat <- scatter3D(M$x, M$y, M$z, pch = "+", cex = 3, colkey = FALSE)
trns3D> # add line
trns3D> XY <- trans3D(x = c(-1, 1), y = c(-1, 1), z = c(-1, 1), pmat = pmat)
trns3D> lines(XY, lwd = 2, col = "blue")
trns3D> ## ========================================================================
trns3D> ## Example 2
trns3D> ## ========================================================================
trns3D>
trns3D> pmat <- perspbox (z = diag(2))
trns3D> XY <- trans3D(x = runif(30), y = runif(30), z = runif(30), pmat = pmat)
trns3D> polygon(XY, col = "darkblue")
plt.pl> # save plotting parameters
plt.pl> pm <- par("mfrow")
plt.pl> pmar <- par("mar")
plt.pl> ## ========================================================================
plt.pl> ## The volcano
plt.pl> ## ========================================================================
plt.pl>
plt.pl> par(mfrow = c(2, 2), mar = c(2, 2, 2, 2))
plt.pl> # The volcano at lower resolution
plt.pl> x <- seq(1, nrow(volcano), by = 2)
plt.pl> y <- seq(1, ncol(volcano), by = 2)
plt.pl> V <- volcano[x,y]
plt.pl> persp3D(z = V)
plt.pl> # rotate
plt.pl> plotdev(theta = 0)
plt.pl> # light and transparence
plt.pl> plotdev(lighting = TRUE, lphi = 90, alpha = 0.6)
plt.pl> # zoom
plt.pl> plotdev(xlim = c(0.2, 0.6), ylim = c(0.2, 0.6), phi = 60)
plt.pl> ## ========================================================================
plt.pl> ## Two spheres
plt.pl> ## ========================================================================
plt.pl>
plt.pl> par(mfrow = c(1, 1), mar = c(0, 0, 0, 0))
plt.pl> # create a sphere
plt.pl> M <- mesh(seq(0, 2*pi, length.out = 30),
plt.pl+ seq(0, pi, length.out = 30))
plt.pl> u <- M$x ; v <- M$y
plt.pl> x <- cos(u)*sin(v)
plt.pl> y <- sin(u)*sin(v)
plt.pl> z <- cos(v)
plt.pl> surf3D(x = 2*x, y = 2*y, z = 2*z,
plt.pl+ colvar = NULL, lighting = TRUE, #plot = FALSE,
plt.pl+ facets = NA, col = "blue", lwd = 5)
plt.pl> surf3D(x, y, z, colvar = NULL, lighting = TRUE,
plt.pl+ col = "red", add = TRUE)
plt.pl> names(getplist())
[1] "type" "plt" "persp" "alpha" "setlim" "xlim"
[7] "ylim" "zlim" "scalefac" "mat" "dot" "imgnr"
[13] "img" "poly"
plt.pl> # plot with different view:
plt.pl> plotdev(phi = 0)
plt.pl> ## Not run:
plt.pl> ##D # will plot same 3-D graph to pdf
plt.pl> ##D pdf(file = "save.pdf")
plt.pl> ##D plotdev()
plt.pl> ##D dev.off()
plt.pl> ## End(Not run)
plt.pl>
plt.pl> ## ========================================================================
plt.pl> ## Two spheres and two planes
plt.pl> ## ========================================================================
plt.pl>
plt.pl> par(mar = c(2, 2, 2, 2))
plt.pl> # equation of a sphere
plt.pl> M <- mesh(seq(0, 2*pi, length.out = 100), -
plt.pl+ seq(0, pi, length.out = 100))
plt.pl> u <- M$x ; v <- M$y
plt.pl> x <- cos(u)*sin(v)
plt.pl> y <- sin(u)*sin(v)
plt.pl> z <- cos(v)
plt.pl> surf3D(x, y, z, colvar = z,
plt.pl+ theta = 45, phi = 20, bty = "b",
plt.pl+ xlim = c(-1.5, 1.5), ylim = c(-1, 2),
plt.pl+ zlim = c(-1.5, 1.5), plot = FALSE)
plt.pl> # add a second sphere, shifted 1 unit to the right on y-axis;
plt.pl> # no facets drawn for this sphere
plt.pl> surf3D (x, y+1, z, colvar = z, add = TRUE,
plt.pl+ facets = FALSE, plot = FALSE)
plt.pl> # define a plane at z = 0
plt.pl> Nx <- 100
plt.pl> Ny <- 100
plt.pl> x <- seq(-1.5, 1.5, length.out = Nx)
plt.pl> y <- seq(-1, 2, length.out = Ny)
plt.pl> image3D (x = x, y = y, z = 0, add = TRUE, colvar = NULL,
plt.pl+ col = "blue", facets = TRUE, plot = FALSE)
plt.pl> # another, small plane at y = 0 - here x and y have to be matrices!
plt.pl> x <- seq(-1., 1., length.out = 50)
plt.pl> z <- seq(-1., 1., length.out = 50)
plt.pl> image3D (x = x, y = 0, z = z, colvar = NULL,
plt.pl+ add = TRUE, col = NA, border = "blue",
plt.pl+ facets = TRUE, plot = TRUE)
plt.pl> ## Not run:
plt.pl> ##D # rotate
plt.pl> ##D for (angle in seq(0, 360, by = 10))
plt.pl> ##D plotdev(theta = angle)
plt.pl> ## End(Not run)
plt.pl>
plt.pl> ## ========================================================================
plt.pl> ## Zooming, rescaling, lighting,...
plt.pl> ## ========================================================================
plt.pl>
plt.pl> par(mfrow = c(2, 2))
plt.pl> # The volcano
plt.pl> x <- seq(1, nrow(volcano), by = 2)
plt.pl> y <- seq(1, ncol(volcano), by = 2)
plt.pl> V <- volcano[x,y]
plt.pl> # plot the volcano
plt.pl> persp3D (x, y, z = V, colvar = V, theta = 10, phi = 20,
plt.pl+ box = FALSE, scale = FALSE, expand = 0.3,
plt.pl+ clim = range(V), plot = FALSE)
plt.pl> # add a plane (image) at z = 170; jetcolored, transparant: only border
plt.pl> image3D(x, y, z = 170, add = TRUE, clim = range(V),
plt.pl+ colvar = V, facets = NA, plot = FALSE, colkey = FALSE)
plt.pl> # add a contour (image) at z = 170; jetcolored,
plt.pl> contour3D(x, y, z = 170, add = TRUE, clim = range(V),
plt.pl+ colvar = V, plot = FALSE, colkey = FALSE)
plt.pl> # plot it -
plt.pl> plot(getplist()) # same as plotdev()
plt.pl> # plot but with different expansion
plt.pl> plotdev(expand = 1)
plt.pl> # other perspective, and shading
plt.pl> plotdev(d = 2, r = 10, shade = 0.3)
plt.pl> # zoom and rotate
plt.pl> plotdev(xlim = c(10, 30), ylim = c(20, 30), phi = 50)
plt.pl> ## ========================================================================
plt.pl> ## Using setplist
plt.pl> ## ========================================================================
plt.pl>
plt.pl> polygon3D(runif(3), runif(3), runif(3))
plt.pl> # retrieve plotting list
plt.pl> plist <- getplist()
plt.pl> names(plist)
[1] "type" "plt" "persp" "alpha" "setlim" "xlim"
[7] "ylim" "zlim" "scalefac" "mat" "dot" "imgnr"
[13] "img" "poly"
plt.pl> plist$poly
$x
[,1]
[1,] 0.3147765
[2,] 0.5926986
[3,] 0.6017956
[4,] NA
$y
[,1]
[1,] 0.60047536
[2,] 0.82227459
[3,] 0.04894024
[4,] NA
$z
[,1]
[1,] 0.6877397
[2,] 0.7195563
[3,] 0.8633404
[4,] NA
$col
[1] "grey"
$border
[1] NA
$lwd
[1] 1
$lty
[1] 1
$alpha
[1] NA
$isimg
[1] 0
$proj
[1] -9.339463
attr(,"class")
[1] "poly"
plt.pl> # change copy of plotting list
plt.pl> plist$poly$col <- "red"
plt.pl> # update internal plotting list
plt.pl> setplist(plist)
plt.pl> # plot updated list
plt.pl> plotdev()
plt.pl> ## ========================================================================
plt.pl> ## Using selectplist
plt.pl> ## ========================================================================
plt.pl>
plt.pl> polygon3D(runif(10), runif(10), runif(10), col = "red",
plt.pl+ alpha = 0.2, plot = FALSE, ticktype = "detailed",
plt.pl+ xlim = c(0,1), ylim = c(0, 1), zlim = c(0, 1))
plt.pl> polygon3D(runif(10)*0.5, runif(10), runif(10), col = "yellow",
plt.pl+ alpha = 0.2, plot = FALSE, add = TRUE)
plt.pl> polygon3D(runif(10)*0.5+0.5, runif(10), runif(10), col = "green",
plt.pl+ alpha = 0.2, plot = FALSE, add = TRUE)
plt.pl> points3D(runif(10), runif(10), runif(10), col = "blue",
plt.pl+ add = TRUE, plot = FALSE)
plt.pl> segments3D(x0 = runif(10), y0 = runif(10), z0 = runif(10),
plt.pl+ x1 = runif(10), y1 = runif(10), z1 = runif(10),
plt.pl+ colvar = 1:10, add = TRUE, lwd = 3)
plt.pl> # retrieve plotting list
plt.pl> plist <- getplist()
plt.pl> # selection function
plt.pl> SS <- function (x, y, z) {
plt.pl+ sel <- rep(TRUE, length.out = length(x))
plt.pl+ sel[x < 0.5] <- FALSE
plt.pl+ return(sel)
plt.pl+ }
plt.pl> # The whole polygon will be removed or kept.
plt.pl> plot(x = selectplist(plist, SS),
plt.pl+ xlim = c(0, 1), ylim = c(0, 1), zlim = c(0, 1))
plt.pl> # restore plotting parameters
plt.pl> par(mfrow = pm)
plt.pl> par(mar = pmar)
ImgOcn> # save plotting parameters
ImgOcn> pm <- par("mfrow")
ImgOcn> mar <- par("mar")
ImgOcn> ## =======================================================================
ImgOcn> ## Images of the hypsometry
ImgOcn> ## =======================================================================
ImgOcn>
ImgOcn> par(mfrow = c(2, 2))
ImgOcn> image2D(Hypsometry, asp = TRUE, xlab = expression(degree*E),
ImgOcn+ ylab = expression(degree*N), contour = TRUE)
ImgOcn> # remove ocean
ImgOcn> zz <- Hypsometry$z
ImgOcn> zz[zz < 0] <- NA
ImgOcn> image2D(zz, x = Hypsometry$x, y = Hypsometry$y, NAcol = "black")
ImgOcn> ## =======================================================================
ImgOcn> ## A short version for plotting the Ocean's bathymetry
ImgOcn> ## =======================================================================
ImgOcn>
ImgOcn> ImageOcean(asp = TRUE, contour = TRUE)
ImgOcn> ImageOcean(col = "white",
ImgOcn+ contour = list(levels = seq(-6000, 0, by = 2000)))
ImgOcn> ## =======================================================================
ImgOcn> ## A complex image of part of the ocean
ImgOcn> ## =======================================================================
ImgOcn>
ImgOcn> # elaborate version
ImgOcn> par(mfrow = c(1, 1), mar = c(2, 2, 2, 2))
ImgOcn> ii <- which(Hypsometry$x > -50 & Hypsometry$x < -20)
ImgOcn> jj <- which(Hypsometry$y > 10 & Hypsometry$y < 40)
ImgOcn> # Draw empty persp box
ImgOcn> zlim <- c(-10000, 0)
ImgOcn> pmat <- perspbox(z = Hypsometry$z[ii, jj],
ImgOcn+ xlab = "longitude", ylab = "latitude", zlab = "depth",
ImgOcn+ expand = 0.5, d = 2, zlim = zlim, phi = 20, theta = 30,
ImgOcn+ colkey = list(side = 1))
ImgOcn> # A function that makes a black panel with grey edge:
ImgOcn> panelfunc <- function(x, y, z) {
ImgOcn+ XY <- trans3D(x, y, z, pmat = pmat)
ImgOcn+ polygon(XY$x, XY$y, col = "black", border = "grey")
ImgOcn+ }
ImgOcn> # left panel
ImgOcn> panelfunc(x = c(0, 0, 0, 0), y = c(0, 0, 1, 1),
ImgOcn+ z = c(zlim[1], zlim[2], zlim[2], zlim[1]))
ImgOcn> # back panel
ImgOcn> panelfunc(x = c(0, 0, 1, 1), y = c(1, 1, 1, 1),
ImgOcn+ z = c(zlim[1], zlim[2], zlim[2], zlim[1]))
ImgOcn> # bottom panel
ImgOcn> panelfunc(x = c(0, 0, 1, 1), y = c(0, 1, 1, 0),
ImgOcn+ z = c(zlim[1], zlim[1], zlim[1], zlim[1]))
ImgOcn> # Actual bathymetry, 2 times increased resolution, with contours
ImgOcn> persp3D(z = Hypsometry$z[ii,jj], add = TRUE, resfac = 2,
ImgOcn+ contour = list(col = "grey", side = c("zmin", "z")),
ImgOcn+ zlim = zlim, clab = "depth, m",
ImgOcn+ colkey = list(side = 1, length = 0.5, dist = -0.1))
ImgOcn> # shorter alternative for same plot, higher resolution
ImgOcn> ## Not run:
ImgOcn> ##D persp3D(z = Hypsometry$z[ii,jj], resfac = 4,
ImgOcn> ##D contour = list(col = "grey", side = c("zmin", "z")),
ImgOcn> ##D zlim = zlim, clab = "depth, m", bty = "bl2",
ImgOcn> ##D xlab = "longitude", ylab = "latitude", zlab = "depth",
ImgOcn> ##D expand = 0.5, d = 2, phi = 20, theta = 30,
ImgOcn> ##D colkey = list(side = 1, length = 0.5, dist = -0.1))
ImgOcn> ## End(Not run)
ImgOcn>
ImgOcn> # reset plotting parameters
ImgOcn> par(mfrow = pm)
ImgOcn> par(mar = mar)
Oxsat> # save plotting parameters
Oxsat> pm <- par("mfrow")
Oxsat> ## ========================================================================
Oxsat> ## plot all surface data
Oxsat> ## ========================================================================
Oxsat>
Oxsat> par(mfrow = c(1, 1))
Oxsat> image2D(z = Oxsat$val[ , , 1], x = Oxsat$lon, y = Oxsat$lat,
Oxsat+ main = "surface oxygen saturation (%) for 2005")
Oxsat> ## ========================================================================
Oxsat> ## plot a selection of latitude-depth profiles; input is an array
Oxsat> ## ========================================================================
Oxsat>
Oxsat> lon <- Oxsat$lon
Oxsat> image2D (z = Oxsat$val, margin = c(2, 3), x = Oxsat$lat,
Oxsat+ y = Oxsat$depth, subset = (lon > 18 & lon < 23),
Oxsat+ ylim = c(5500, 0), NAcol = "black", zlim = c(0, 110),
Oxsat+ xlab = "latitude", ylab = "depth, m")
Oxsat> ImageOcean()
Oxsat> abline ( v = lon[lon > 18 & lon < 23])
Oxsat> ## ========================================================================
Oxsat> ## plot with slices
Oxsat> ## ========================================================================
Oxsat>
Oxsat> par(mfrow = c(1, 1))
Oxsat> ii <- which (Oxsat$lon > -90 & Oxsat$lon < 90)
Oxsat> jj <- which (Oxsat$lat > 0 & Oxsat$lat < 90)
Oxsat> xs <- Oxsat$lon[ii[length(ii)]] # E boundary
Oxsat> ys <- Oxsat$lat[jj[1]] # S boundary
Oxsat> slice3D(colvar = Oxsat$val[ii,jj,], x = Oxsat$lon[ii],
Oxsat+ y = Oxsat$lat[jj], z = -Oxsat$depth,
Oxsat+ NAcol = "black", xs = xs, ys = ys, zs = 0,
Oxsat+ theta = 35, phi = 50, colkey = list(length = 0.5),
Oxsat+ expand = 0.5, ticktype = "detailed",
Oxsat+ clab = "%", main = "Oxygen saturation",
Oxsat+ xlab = "longitude", ylab = "latitude", zlab = "depth")
Oxsat> # restore plotting parameters
Oxsat> par(mfrow = pm)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.