// wrapper to GSL spline
// spline type "cspline" is nat. cubic spline
// spline type "akima" or "Akima" is Akima spline
// C-style string comparison!
#ifndef GSL_SPLINE_HPP
#define GSL_SPLINE_HPP
#include <iostream>
#include <cstring>
#include <vector>
#include <gsl/gsl_errno.h>
#include <gsl/gsl_spline.h>

using namespace std;

void spline(const vector<double>& x,const vector<double>& y,
	    vector<double>& xx, vector<double>& yy,const char* type){  
  gsl_interp_accel* acc  = gsl_interp_accel_alloc();
  gsl_spline* spline;
  int choose=0;
  if(strcmp(type,"akima")==0) choose=1;
  if(strcmp(type,"Akima")==0) choose=1;
  switch (choose){
  case 0:
    spline = gsl_spline_alloc(gsl_interp_cspline, x.size());
    break;
  case 1: 
    spline = gsl_spline_alloc(gsl_interp_akima, x.size());
    break;
  }
  gsl_spline_init(spline, &x[0], &y[0], x.size());
  vector<double>::iterator posx,posy;
  posy = yy.begin();
  for(posx=xx.begin();posx!=xx.end();++posx){
    *posy++= gsl_spline_eval(spline, *posx, acc);
  }
  gsl_spline_free(spline);
  gsl_interp_accel_free(acc);
}
#endif