Commit 4fc50385 authored by agebhard's avatar agebhard
Browse files

final changes before upload of akima 0.4

parent 8c525f1e
Package: akima
Version: 0.4-1
Version: 0.4-3
Title: Interpolation of irregularly spaced data
Author: Fortran code by H. Akima
R port by Albrecht Gebhardt <albrecht.gebhardt@uni-klu.ac.at>
......
aspline <- function(x, y, xout, n = 50, ties=mean, algorithm="original", np=3) {
aspline <- function(x, y=NULL, xout, n = 50, ties=mean, method="original", degree=3) {
if (! algorithm %in% c("original", "improved")) stop(paste("unknown algorithm:", algorithm))
if (! method %in% c("original", "improved")) stop(paste("unknown method:", method))
x <- xy.coords(x, y) # -> (x,y) numeric of same length
y <- x$y
......@@ -29,24 +29,38 @@ aspline <- function(x, y, xout, n = 50, ties=mean, algorithm="original", np=3) {
stop("need at least two non-NA values to interpolate")
if (missing(xout)) {
if (n <= 0)
stop("approx requires n >= 1")
stop("aspline requires n >= 1")
xout <- seq(x[1], x[nx], length = n)
}
nout <- length(xout)
yout <- numeric(nout)
if (algorithm == "improved") {
y <- .Fortran("uvip3p", as.integer(np),
err <- 0
if (method == "improved") {
ret <- .Fortran("uvip3p", as.integer(degree),
as.integer(nx), as.double(x), as.double(y),
as.integer(nout), as.double(xout), yout=as.double(yout),
PACKAGE="akima")$yout
as.integer(nout), as.double(xout), yout=as.double(yout),
err=as.integer(err), PACKAGE="akima")
} else {
y <- .Fortran("intrpl",
ret <- .Fortran("intrpl",
as.integer(nx), as.double(x), as.double(y),
as.integer(nout), as.double(xout), yout=as.double(yout),
PACKAGE="akima")$yout
as.integer(nout), as.double(xout), yout=as.double(yout),
err=as.integer(err), PACKAGE="akima")
}
list(x = xout, y = y)
err <- max(0, min(ret$err, 10)) # 10 is maximum error code
if (err > 0) {
## if the error handling befor .Fortran is correct
## the following lines should never be called
errl <- c("Insufficient number of data points",
"No desired points",
"Two data points identical or out of sequence",
"undefined error code",
"undefined error code",
"identical x values",
"x values out of sequence",
"undefined error code",
"undefined error code",
"internal error in Akima Fortran code")
warning(errl[err])
}
list(x = xout, y = ret$yout)
}
\name{aspline}
\alias{aspline}
%- Also NEED an '\alias' for EACH other topic documented here.
\title{Univariate Akima interpolation}
\description{
The function returns a list of points which smoothly
interpolate given data points, similar to a curve drawn by hand.}
\usage{
aspline(x, y, xout, n = 50, ties = mean, algorithm="original", np=3) }
%- maybe also 'usage' for other objects documented here.
aspline(x, y=NULL, xout, n = 50, ties = mean, method="original", degree=3) }
\arguments{
\item{x, y}{vectors giving the coordinates of the points to be
......@@ -23,9 +21,9 @@
\item{ties}{Handling of tied \code{x} values. Either a function
with a single vector argument returning a single number result or
the string \code{"ordered"}.}
\item{algorithm}{either \code{"original"} method after Akima (1970) or
\item{method}{either \code{"original"} method after Akima (1970) or
\code{"improved"} method after Akima (1991)}
\item{np}{im improved algorithm is selected: degree of the
\item{degree}{if improved algorithm is selected: degree of the
polynomials for the interpolating function}
}
......@@ -46,18 +44,14 @@
\references{
Akima, H. (1970) A new method of interpolation
and smooth curve fitting based on local procedures,
J. ACM 17(4), 589-602
J. ACM \bold{17}(4), 589-602
Akima, H. (1991) A Method of Univariate Interpolation that Has
the Accuracy of a Third-degree Polynomial. ACM Transactions on
Mathematical Software, 17(3), 341-366.
Mathematical Software, \bold{17}(3), 341-366.
}
\author{Thomas Petzoldt}
%%\note{ ~~further notes~~ }
%% ~Make other sections like Warning with \section{Warning }{....} ~
\seealso{\code{\link{approx}}, \code{\link{spline}}}
\examples{
......@@ -70,8 +64,8 @@ plot(x, y, ylim=c(-3, 3), xlim=range(xnew))
lines(spline(x, y, xmin=min(xnew), xmax=max(xnew), n=200), col="blue")
lines(aspline(x, y, xnew), col="red")
lines(aspline(x, y, xnew, algorithm="improved"), col="black", lty="dotted")
lines(aspline(x, y, xnew, algorithm="improved", np=10), col="green", lty="dashed")
lines(aspline(x, y, xnew, method="improved"), col="black", lty="dotted")
lines(aspline(x, y, xnew, method="improved", degree=10), col="green", lty="dashed")
## irregular spaced data
x <- sort(runif(10, max=10))
......@@ -82,8 +76,8 @@ plot(x, y, ylim=c(-3, 3), xlim=range(xnew))
lines(spline(x, y, xmin=min(xnew), xmax=max(xnew), n=200), col="blue")
lines(aspline(x, y, xnew), col="red")
lines(aspline(x, y, xnew, algorithm="improved"), col="black", lty="dotted")
lines(aspline(x, y, xnew, algorithm="improved", np=10), col="green", lty="dashed")
lines(aspline(x, y, xnew, method="improved"), col="black", lty="dotted")
lines(aspline(x, y, xnew, method="improved", degree=10), col="green", lty="dashed")
## an example of Akima, 1991
x <- c(-3, -2, -1, 0, 1, 2, 2.5, 3)
......@@ -93,8 +87,8 @@ plot(x, y, ylim=c(-3, 3))
lines(spline(x, y, n=200), col="blue")
lines(aspline(x, y, n=200), col="red")
lines(aspline(x, y, n=200, algorithm="improved"), col="black", lty="dotted")
lines(aspline(x, y, n=200, algorithm="improved", np=10), col="green", lty="dashed")
lines(aspline(x, y, n=200, method="improved"), col="black", lty="dotted")
lines(aspline(x, y, n=200, method="improved", degree=10), col="green", lty="dashed")
}
\keyword{arith}
\keyword{dplot}
\ No newline at end of file
\keyword{dplot}
......@@ -3,7 +3,7 @@ C
C ALGORITHM 433 COLLECTED ALGORITHMS FROM ACM.
C ALGORITHM APPEARED IN COMM. ACM, VOL. 15, NO. 10,
C P. 914.
SUBROUTINE INTRPL(L,X,Y,N,U,V)
SUBROUTINE INTRPL(L,X,Y,N,U,V,ERR)
C INTERPOLATION OF A SINGLE-VALUED FUNCTION
C THIS SUBROUTINE INTERPOLATES, FROM VALUES OF THE FUNCTION
C GIVEN AS ORDINATES OF INPUT DATA POINTS IN AN X-Y PLANE
......@@ -25,8 +25,9 @@ C (ABSCISSAS) OF DESIRED POINTS
C THE OUTPUT PARAMETER IS
C V = ARRAY OF DIMENSION N WHERE THE INTERPOLATED Y
C VALUES (ORDINATES) ARE TO BE DISPLAYED
C ERR = ERROR CODE (added by ThPe)
C DECLARATION STATEMENTS
INTEGER L, N
INTEGER L, N, ERR
DOUBLE PRECISION X(L),Y(L),U(N),V(N)
EQUIVALENCE (P0,X3),(Q0,Y3),(Q1,T3)
DOUBLE PRECISION M1,M2,M3,M4,M5
......@@ -34,7 +35,7 @@ C DECLARATION STATEMENTS
EQUIVALENCE (UK,DX),(IMN,X2,A1,M1),(IMX,X5,A5,M5),
1 (J,SW,SA),(Y2,W2,W4,Q2),(Y5,W3,Q3)
C PRELIMINARY PROCESSING
10 L0=L
10 L0=L
LM1=L0-1
LM2=LM1-1
LP1=L0+1
......@@ -154,22 +155,14 @@ C COMPUTATION OF THE POLYNOMIAL
80 V(K)=Q0+DX*(Q1+DX*(Q2+DX*Q3))
RETURN
C ERROR EXIT
90 WRITE (2090)
GO TO 99
91 WRITE (2091)
GO TO 99
95 WRITE (2095)
GO TO 97
96 WRITE (2096)
97 WRITE (2097) I,X(I)
99 WRITE (2099) L0,N0
90 ERR = 1
RETURN
91 ERR = 2
RETURN
95 ERR = 6
RETURN
96 ERR = 7
RETURN
99 ERR = 10
RETURN
C FORMAT STATEMENTS
2090 FORMAT(1X/22H *** L = 1 OR LESS./)
2091 FORMAT(1X/22H *** N = 0 OR LESS./)
2095 FORMAT(1X/27H *** IDENTICAL X VALUES./)
2096 FORMAT(1X/33H *** X VALUES OUT OF SEQUENCE./)
2097 FORMAT(6H I =,I7,10X,6HX(I) =,E12.3)
2099 FORMAT(6H L =,I7,10X,3HN =,I7/
1 36H ERROR DETECTED IN ROUTINE INTRPL)
END
......@@ -4,7 +4,7 @@ C ALGORITHM 697 , COLLECTED ALGORITHMS FROM ACM.
C THIS WORK PUBLISHED IN TRANSACTIONS ON MATHEMATICAL SOFTWARE,
C VOL. 17, NO. 3, SEPTEMBER, 1991, PP. 367.
C
SUBROUTINE UVIP3P(NP,ND,XD,YD,NI,XI, YI)
SUBROUTINE UVIP3P(NP,ND,XD,YD,NI,XI, YI, ERR)
C
C Univariate Interpolation (Improved Akima Method)
C
......@@ -48,6 +48,7 @@ C
C The output argument is
C YI = array of dimension NI, where the ordinates of the
C desired points are to be stored.
C ERR = ERROR CODE (added by ThPe)
C
C If an integer value smaller than 3 is given to the NP argument,
C this subroutine assumes NP = 3.
......@@ -388,20 +389,12 @@ C (Cubic interpolation and linear extrapolation)
C End of Special Case 3
RETURN
C Error exit
90 WRITE (*,99090) ND
GO TO 99
91 WRITE (*,99091) NI
GO TO 99
92 WRITE (*,99092) ID,XD(ID-1),XD(ID)
99 WRITE (*,99099)
90 ERR=1
RETURN
91 ERR=2
RETURN
92 ERR=3
RETURN
99 ERR=10
RETURN
C Format statements for error messages
99090 FORMAT (1X/ ' *** Insufficient data points.'/
1 7X,'ND =',I3)
99091 FORMAT (1X/ ' *** No desired points.'/
1 7X,'NI =',I3)
99092 FORMAT (1X/ ' *** Two data points identical or out of ',
1 'sequence.'/
2 7X,'ID, XD(ID-1), XD(ID) =',I5,2F10.3)
99099 FORMAT (' Error detected in the UVIP3P subroutine'/)
END
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment