bicubic_interp.java

//u and v are the two parameters of the surface between 0 and numPtsX
//ignoring boundary conditions
PVector bicubic(PVector[][] pts, float u, float v) {
  //get indices
  int uInd1 = int(u);
  int uInd0 = uInd1-1;
  int uInd2 = uInd1+1;
  int uInd3 = uInd1+2;
  float uVal = u-uInd1; //between 0..1

  int vInd1 = int(v);
  int vInd0 = vInd1-1;
  int vInd2 = vInd1+1;
  int vInd3 = vInd1+2;
  float vVal = v-vInd1; //between 0..1
  
  //do first four interpolations
  PVector p0 = cubic(new PVector[]
     {pts[uInd0][vInd0],pts[uInd1][vInd0],pts[uInd2][vInd0],pts[uInd3][vInd0]}
     ,uVal);
  PVector p1 = cubic(new PVector[]
     {pts[uInd0][vInd1],pts[uInd1][vInd1],pts[uInd2][vInd1],pts[uInd3][vInd1]}
     ,uVal);
  PVector p2 = cubic(new PVector[]
     {pts[uInd0][vInd2],pts[uInd1][vInd2],pts[uInd2][vInd2],pts[uInd3][vInd2]}
     ,uVal);
  PVector p3 = cubic(new PVector[]
     {pts[uInd0][vInd3],pts[uInd1][vInd3],pts[uInd2][vInd3],pts[uInd3][vInd3]}
     ,uVal);

  return cubic(new PVector[]{p0,p1,p2,p3},vVal);
}

PVector cubic(PVector[] pts, float t) {
  PVector pt = new PVector();
  pt.add(PVector.mult(pts[0],-.5*t+t*t+.5*t*t*t));
  pt.add(PVector.mult(pts[1],1-2.5*t*t+1.5*t*t*t));
  pt.add(PVector.mult(pts[2],.5*t+2*t*t-1.5*t*t*t));
  pt.add(PVector.mult(pts[3],-.5*t*t+.5*t*t*t));
}