Listing 3: Piecewise linear interpolation

#include "interpolate.h"

/* returns piecewise-linear approximation of f(x) at X
double pwlinear
   (
   const double *x,  /* array of strictly increasing points */
   const double *fx, /* array of known values of f(x)       */
   int n,            /* number of values of x and fx        */
   double X          /* value of x at which to compute f(x) */
   ) 
{

   int                k;

   /*
   ** Locate the interval on which X lies. If this interval is 
   ** not found by k = n - 3, [x[n - 2], x[n - 1]] is used.
   */

   for (k = 0; k < n - 2; k++) {

      if (X < x[k + 1])
         break;

   }

   /*
   ** Compute p1(x), the linear interpolating polynomial, 
   ** to obtain a value for f(x) at X.      
   */

   return fx[k] + 
             (((fx[k+1] - fx[k])/(x[k+1] - x[k]))*(X - x[k]));

}