Many thanks to Martin Maechler and Jim Garrett.
They proposed almost the same algorithm to solve the problem (how to draw an
interpolating curve through a set of points).
Martin Maechler wrote a nice function. A particular problem arises if the shape
must be a closed shape.
So I added some elements to obtain two functions, one for an open shape, one for
a closed. They can be used as an alternative of lines() and polygon(): sometimes
the angles can be seen as artificial, it is a matter of visual semantics.
These functions are not very smart, many improvements are surely possible.
****************************************************************************************************
spllines <- function(x, y = NULL, n = 20, type="l",
col=par("col"), lty=par("lty"), ... ){
## Purpose: Interpolating spline *curve* (not function)
## for an open shape (missing values not accepted)
## ---------------------------------------------------------------
## Arguments: (almost) the same as for lines()
## ---------------------------------------------------------------
## Author: Martin Maechler, Date: 6 Dec 2002, 12:34
xy <- cbind(x,y)
nP <- length(xy[,1])
if((nP < 3) || (is.na(any(xy[,1]))))
return(list(x=numeric(0), y=numeric(0)))
## else :
nP <- length(xy[,1])
z <- n*(nP-1)
i <- 1:nP
Sx <- splinefun(i, xy[,1], method = "natural")
Sy <- splinefun(i, xy[,2], method = "natural")
ti <- seq(1, nP, length = z)
opspl <- cbind(x = Sx(ti), y = Sy(ti))
lines(opspl, type=type, col=col, lty=lty, ... )
}
****************************************************************************
splpolygon <- function(x, y = NULL, n = 20, density=NULL, angle=45,
border=NULL, col=NA, lty=NULL, xpd=NULL, ... ){
## Purpose: Interpolating spline *curve* (not function)
## for a closed shape (missing values not accepted)
## ---------------------------------------------------------------
## Arguments: (almost) the same as for polygon()
## ---------------------------------------------------------------
## Author: Martin Maechler, Date: 6 Dec 2002, 12:34
xy <- cbind(x,y)
nP <- length(xy[,1])
if((nP < 3) || (is.na(any(xy[,1]))))
return(list(x=numeric(0), y=numeric(0)))
## else :
if(xy[1,] == xy[nP,]) xy <- rbind(xy[(nP-1),],xy[nP,],xy,xy[1,],xy[2,])
else xy <- rbind(xy[(nP-2),],xy[(nP-1),],xy[nP,],xy,xy[1,],xy[2,])
nP <- length(xy[,1])
z <- n*(nP-1)
i <- 1:nP
Sx <- splinefun(i, xy[,1], method = "natural")
Sy <- splinefun(i, xy[,2], method = "natural")
ti <- seq(1, nP, length = z)
clspl <- cbind(x = Sx(ti), y = Sy(ti))
clspl <- clspl[((n*2)+1):(z-(n*2)),]
polygon(clspl, density=density, angle=angle, border=border, col=col, lty=lty,
xpd=xpd, ...)
}
**************************************************************
Alain GUERREAU CNRS-Paris
[[alternate HTML version deleted]]