This library contains an R implementation of the S-Plus function interp().
I used Fortran code for cubic spline interpolation from H. Akima available
at http://www.netlib.org/toms .
There are two version, one from 1978 (AMS code 526, used in SPlus 3.x) and a
newer one (AMS 761) from 1996. The original versions are in the directory orig/.
I only had to write a subroutine idptli.f for linear interpolation within
triangles from the delauney triangulation.
I added a function interpp(), which evaluates the interpolated function at
arbitraryly choosen points and generates no regular grid as interp() does.
interp() actually calls interp.new() or interp.old() depending on its arguments,
linear interpolation is done by interp.old(), spline interpolation by interp.new
.
UPDATE (13.12.2004):
Thomas Petzoldt provided an R interface to
Akimas univariate spline interpolation functions (aspline) which has been
included into this library.
The original papers from H. Akima are available (access is not free, depends on your
ACM subscription):
TOMS 433: http://dl.acm.org/citation.cfm?id=355604.355605
TOMS 526: http://dl.acm.org/citation.cfm?id=355780.355787
TOMS 697: http://dl.acm.org/citation.cfm?id=114697.116808
TOMS 761: http://dl.acm.org/citation.cfm?id=232826.232856
Original Fortran code is freely available at Netlib/TOMS -- Collected Algorithms
of the ACM at:
http://www.netlib.org/toms/433
http://www.netlib.org/toms/526
http://www.netlib.org/toms/697
http://www.netlib.org/toms/761
------------------------------------------------------------------------------
Albrecht Gebhardt email: albrecht.gebhardt@aau.at
Institut fuer Mathematik Tel. : (++43 463) 2700/3118
Universitaet Klagenfurt Fax : (++43 463) 2700/3198
Universitaetsstr. 65-67 WWW : http://www.stat.aau.at/~agebhard
A-9020 Klagenfurt, Austria
------------------------------------------------------------------------------
Details:
See ChangeLog for details about further versions
Version 0.1-1
I used the S-Plus version of interp() (S-Plus 3.3 f. Win 3.11) as starting
point. Then I searched for Akimas idsfft function, which is referenced in
interp().
I found the appropriate Fortran code (1978) in the ACM Collected
Algorithms archive under
http://www.netlib.org/toms/526
(also included here) and splitted it into the Fortran files under src/ .
The test driver ttidbs was also included in older versions and could be
compiled with "make test", now it is removed here.
However, it seems that this code differs a little bit from that used in S-Plus:
It implements "only" bicubic spline interpolation, no linear interpolation as
in S-Plus.
So I modified IDSFFT and added a subroutine IDPTLI which does linear
interpolation within the triangles, generated by IDTANG.
Further changes are
+ REAL -> DOUBLE PRECISION
+ static DIMENSIONs replaced with dynamic
+ option to toggle extrapolation outside of convex hull added in IDSFFT and
IDBVIP.
Because I don't know how to generate NAs in Fortran (I use g77 on Linux),
I added a logical array MISSI, that indicates if a returned value should
be NA. These values will be set to NA after the Fortran call.
+ option to handle duplicate data points added (according to an example in the
S-Plus help page)
+ man pages converted and rewritten
+ data set akima (from S-Plus) added
+ function interpp() added, it evaluates the interpolated function at
arbitraryly choosen points and generates no regular grid as interp() does.
There where some problems with interpp() when using the Fortran version:
- it crashes when compiled with "g77 -O2 -fpic" and called with more than
one output point.
- compilation with "g77 -g -fpic" fails (see src/Makefile for details)
- compilation with "g77 -g" works (and no crashes occur)
These problems do not occur in the C Version (generated by f2c), so it would
be better to use only the src-c tree.
After I finished the above steps I found a more recent version (ACM 761, 1996)
of Akimas interpolation code which uses the tripack package (also available at
ACM as algorithm no. 751) for triangulation and now I'm trying to use it for
the next version of interp().