/*    ATK IV - Numerical Programming (2006)      */
/*  Excercise 5.2b -  Cubic Spline Interpolation */

#include<stdio.h>
#include<math.h>
#define NRANSI
#include<nr.h>
#define Pi 3.141592
#define N 20

float f(float x){
  if(x<=-1.0) return x+3.0;
  else if(x<=1) return -2.0*x;
  else return x-3.0;
}

int main(void){
  float xi[N+1],yi[N+1],y2[N+1];
  int i;
  float step;
  float a=-4.0;
  float b=4.0;
  float z=a;
  float y;

  step = (b-a)/(N-1);

  for(i=1;i<=N;i++){
    xi[i]=z;
    yi[i]=f(z);
    z+=step;
  }
  
  spline(xi,yi,N,1,1,y2);

  for(i=1;i<4*N;i++){
    splint(xi,yi,y2,N,a+(i-1)*step/4.0,&y);
    printf("%f %f\n",y,f(a+(i-1)*step/4.0));
    }

  return 0;
}