## -----------------------------------------------------------------------------
z <- function(x, y) { -12 * x^3 + 10 * x^2 * y - 14 * x * y^2 + 25 * y^3 + 50 }
z(.5, 0.8)


## -----------------------------------------------------------------------------
zf <- function(model, xy) {
   x <- xy[,1]
   y <- xy[,2]
   z(x, y)
}


## ----fields5------------------------------------------------------------------
library(terra)
r <- rast(xmin=0.5, xmax=1.4, ymin=0.6, ymax=1.5, ncol=9, nrow=9, crs="")
z <- interpolate(r, model=NULL, fun=zf)
names(z) <- 'z'

vt <- persp(z, theta=30, phi=30, ticktype='detailed', expand=.8)


## ----fields6------------------------------------------------------------------
pts <- as.data.frame(z, xy=TRUE)
pt <- trans3d(pts[,1], pts[,2], pts[,3], vt)
plot(pt, col=rainbow(9, .75, start=.2)[round(pts[,3]/10)-2], pch=20, cex=2)


## ----fields10-----------------------------------------------------------------
library(rspat)
d <- spat_data("precipitation")
head(d)
d$prec <- rowSums(d[, c(6:17)])
plot(sort(d$prec), ylab='Annual precipitation (mm)', las=1, xlab='Stations')


## ----fields15-----------------------------------------------------------------
dsp <- vect(d, c("LONG", "LAT"), crs= "+proj=longlat +datum=WGS84")
CA <- spat_data("counties")

cuts <- c(0,200,300,500,1000,3000)
blues <- colorRampPalette(c('yellow', 'orange', 'blue', 'dark blue'))
plot(CA, col="light gray", lwd=4, border="white")
plot(dsp, "prec", breaks=cuts, col=blues(5), pch=20, cex=1.5, add=TRUE, plg=list(x="topright"))
lines(CA, lwd=1.5)