From 83119fe9bb0c9c2a70ea5397aca9968bde7ffa07 Mon Sep 17 00:00:00 2001 From: tlatorre Date: Fri, 31 Aug 2018 10:33:24 -0500 Subject: add interp1d function to do fast interpolation when the x values are evenly spaced --- src/test.c | 47 ++++++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 46 insertions(+), 1 deletion(-) (limited to 'src/test.c') diff --git a/src/test.c b/src/test.c index cc054c0..2ea73d1 100644 --- a/src/test.c +++ b/src/test.c @@ -11,6 +11,7 @@ #include "path.h" #include "random.h" #include "pdg.h" +#include typedef int testFunction(char *err); @@ -396,6 +397,49 @@ int test_path(char *err) return 0; } +int test_interp1d(char *err) +{ + /* Tests the interp1d() function. */ + size_t i; + double xp[100], yp[100], range, x, y, expected; + gsl_interp_accel *acc; + gsl_spline *spline; + + range = 1.0; + + init_genrand(0); + + for (i = 0; i < sizeof(xp)/sizeof(xp[0]); i++) { + xp[i] = range*i/(100-1); + yp[i] = genrand_real2(); + } + + acc = gsl_interp_accel_alloc(); + spline = gsl_spline_alloc(gsl_interp_linear,100); + gsl_spline_init(spline,xp,yp,100); + + for (i = 0; i < 1000; i++) { + x = genrand_real2()*range; + expected = gsl_spline_eval(spline,x,acc); + y = interp1d(x,xp,yp,100); + if (!isclose(y, expected, 1e-9, 1e-9)) { + sprintf(err, "interp1d(%.2g) returned %.5g, but expected %.5g", x, y, expected); + goto err; + } + } + + gsl_interp_accel_free(acc); + gsl_spline_free(spline); + + return 0; + +err: + gsl_interp_accel_free(acc); + gsl_spline_free(spline); + + return 1; +} + struct tests { testFunction *test; char *name; @@ -408,7 +452,8 @@ struct tests { {test_muon_get_T, "test_muon_get_T"}, {test_logsumexp, "test_logsumexp"}, {test_sno_charge_integral, "test_sno_charge_integral"}, - {test_path, "test_path"} + {test_path, "test_path"}, + {test_interp1d, "test_interp1d"} }; int main(int argc, char **argv) -- cgit