This package provides a minimal set of interpolation methods
(piecewise constant, linear and spline) designed to be compatible with
R’s approx and spline functions but callable
from C. It will only be of interest to people writing C or C++
packages.
The package is designed to be used with R’s LinkingTo:
support and is header only. This is a somewhat awkward situation for C
(rather than C++). The approach is the same as taken by ring.
In your DESCRIPTION, add a line
LinkingTo: cinterpolate. You also need to add
cinterpolate to Imports and ensure that the
package is loaded before use, as it uses R’s R_GetCCallable
interface to call the actual interpolation functions.
In your src/ directory, add a file
cinterpolate.c containing just the line
#include <cinterpolate/cinterpolate.c>.
Anywhere in your C (or C++) code you want to use the
interpolation code buffer, include the line
#include <cinterpolate/cinterpolate.h> to include
the prototypes and use the interface as described below.
(I am not sure what the best practice way of doing this with a
standalone shared library compiled with R CMD SHLIB is
though; probably best to make a package.)
There are only three functions in the cinterpolate API;
one to build the object (cinterpolate_alloc), one for
carrying out interpolation (cinterpolate_eval) and one for
freeing the object after calculations have been run
(cinterpolate_free). If you allocate a cinterpolate object
then you are responsible for freeing it (even on error elsewhere in
code). Not doing this will cause leaks.
#ifndef CINTERPOLTE_CINTERPOLATE_H_
#define CINTERPOLTE_CINTERPOLATE_H_
#include <stddef.h> // size_t
#include <stdbool.h> // bool
// Allow use from C++
#ifdef __cplusplus
extern "C" {
#endif
// There are only three functions in the interface; allocation,
// evaluation and freeing.
// Allocate an interpolation object.
//
// type: The mode of interpolation. Must be one of "constant",
// "linear" or "spline" (an R error is thrown if a different
// value is given).
//
// n: The number of `x` points to interpolate over
//
// ny: the number of `y` points per `x` point. This is 1 in the
// case of zimple interpolation as used by Rs `interpolate()`
//
// x: an array of `x` values of length `n`
//
// y: an array of `ny` sets of `y` values. This is in R's matrix
// order (i.e., the first `n` values are the first series to
// interpolate over).
//
// fail_on_extrapolate: if true, when an extrapolation occurs throw
// an error; if false return NA_REAL
//
// auto_free: automatically clean up the interpolation object on
// return to R. This uses `R_alloc` for allocations rather than
// `Calloc` so freeing will always happen (even on error
// elsewhere in the code). However, this prevents returning back
// a pointer to R that will last longer than the call into C
// code.
//
// The return value is an opaque pointer that can be passed through to
// `cinterpolate_eval` and `cinterpolate_free`
void *cinterpolate_alloc(const char *type, size_t n, size_t ny,
double *x, double *y, bool fail_on_extrapolate,
bool auto_free);
// Evaluate the interpolated function at a new `x` point.
//
// x: A new, single, `x` point to interpolate `y` values to
//
// obj: The interpolation object, as returned by `cinterpolate_alloc`
//
// y: An array of length `ny` to store the interpolated values
//
// The return value is 0 if the interpolation is successful (with x
// lying within the range of values that the interpolation function
// supports), -1 otherwise
int cinterpolate_eval(double x, void *obj, double *y);
// Clean up all allocated memory
//
// obj: The interpolation object, as returned by `cinterpolate_alloc`
void cinterpolate_free(void *obj);
#ifdef __cplusplus
}
#endif
#endifA complete example of use is included in the package as
system.file("example", package = "cinterpolate").
The DESCRIPTION looks like
Package: example
Title: Testing 'cinterpolate' Package Use
Version: 0.0.1
Description: Testing that using 'cinterpolate' from another package works.
License: CC0
Authors@R: person("Rich", "FitzJohn", role = c("aut", "cre"),
email = "rich.fitzjohn@gmail.com")
LinkingTo: cinterpolate
Imports: cinterpolate
Suggests: testthatNote the use of LinkingTo: and Imports:
here.
The NAMESPACE file ensures that the package’s shared
library is loaded (useDynLib(example)) and that
cinterpolate’s functions will be available by importing the
package import(cinterpolate)
(importFrom(cinterpolate, interpolate_function) would also
be fine).
useDynLib(example)
import(cinterpolate)The actual usage from C looks like:
#include <cinterpolate/cinterpolate.c>
#include <stdbool.h>
#include <R.h>
#include <Rinternals.h>
#include <R_ext/Utils.h>
SEXP test(SEXP r_x, SEXP r_y, SEXP r_xout, SEXP r_type) {
const char * type = CHAR(STRING_ELT(r_type, 0));
bool is_matrix = isMatrix(r_y);
size_t n = length(r_x);
size_t ny = is_matrix ? ncols(r_y) : 1;
void *obj = cinterpolate_alloc(type, n, ny, REAL(r_x), REAL(r_y),
false, true);
size_t nxout = length(r_xout);
SEXP r_yout = PROTECT(is_matrix ?
allocMatrix(REALSXP, ny, nxout) :
allocVector(REALSXP, nxout));
double *xout = REAL(r_xout), *yout = REAL(r_yout);
for (size_t i = 0; i < nxout; ++i) {
cinterpolate_eval(xout[i], obj, yout);
yout += ny;
}
UNPROTECT(1);
return r_yout;
}Because this is the only use of cinterpolate in the package, we
can directly include cinterpolate/cinterpolate.c
The interpolate_alloc function is used with
auto_clean = true so there is no use of
interpolate_free - because this interpolator only needs to
survive for this single C function this method of cleanup is probably
better
There is a single allocation of an interpolation object but several calls to interpolate new values.