C implementation of algorithm 2.1 and 2.2 from *The Nurbs Book*

SplineBasis.cpp

#include <stdio.h>
#include <stdlib.h>

int findSpan( int n, int p, double u, double U[]) {
  if (u == U[n+1]) return(n);
  int pos = n/2;
  int step = (pos < 4) ? 1 : pos / 2;
  while (u < U[pos] || u >= U[pos+1]) {
    pos += step * (( u >= U[pos]) * 2 - 1);
    step = (step + 1)/2;
  }
  return pos;
}

void baseFuncs( int i, double u, int p, double U[], double N[]) {
  double left[p+1];
  double right[p+1];
  N[0]=1.0;
  for (int j=1; j<=p; ++j) {
    left[j]  = u - U[i+1-j];
    right[j] = U[i+j] - u;
    double saved = 0.0;
    for (int r=0; r<j; ++r) {
      double temp= N[r]/(right[r+1]+left[j-r]);
      N[r] = saved+right[r+1]*temp;
      saved = left[j-r]*temp;
    }
    N[j] = saved;
  }
}

int main() {
  double U[] = {0, 0, 0.2, 0.3, 0.6, 0.9, 1, 1};
  const int p = 2;
  const int sU = sizeof U / sizeof U[0];
  const int n = sU - p - 1;
  double N[p+1];

//  FILE * outf = stdout; 
  FILE * outf = fopen("base.txt","w");
  for ( double u = 0; u < 1; u += 1e-3 ) {
    fprintf ( outf, "%.6f", u );
    int i = findSpan(n,p,u,U);
    baseFuncs( i, u, p, U, N);
    for (int x = p-1; x < i; ++x) fprintf ( outf, " %.3f", 0. );
    for (int x = 0; x < p+1; ++x) fprintf ( outf, " %.3f", N[x] );
    for (int x = i + p + 1; x < sU; ++x) fprintf ( outf, " %.3f", 0. );
    fprintf ( outf, "\n" );
  }

  if (outf && outf != stdout) fclose(outf);
  return 0;
}