class BezierPatch{
  private float[][][] patch;
  private float uStep=0.125, vStep=0.125;
  //--------------
  BezierPatch(){
    patch = new float[][][]{
    new float[][]{ {0, 0,0},{30, 0,0},{60, 0,0},{90, 0,0} },
    new float[][]{ {0,30,0},{30,30,0},{60,30,0},{90,30,0} },
    new float[][]{ {0,60,0},{30,60,0},{60,60,0},{90,60,0} },
    new float[][]{ {0,90,0},{30,90,0},{60,90,0},{90,90,0} }
    };
  }
  //--------------
  // supply an array with 16 points in 3D space, each represented by a float array of length 3
  BezierPatch(float[][] p){
    this( new float[][][]{
    new float[][]{ p[0],p[1],p[2],p[3] },
    new float[][]{ p[4],p[5],p[6],p[7] },
    new float[][]{ p[8],p[9],p[10],p[11] },
    new float[][]{ p[12],p[13],p[14],p[15] }
    } );
  }
  //--------------
  // supply an 4*4 array with 16 points in 3D space, each represented by a float array of length 3
  BezierPatch(float[][][] p){
    patch = new float[4][4][3];
    for(int i=0 ; i<4 ; i++ )
      for(int j=0 ; j<4 ; j++ )
        for(int c=0 ; c<3 ; c++ )
          patch[i][j][c] = p[i][j][c];
  }
  //--------------
  void setDetail(int res){
    setDetail(res,res);
  }
  //--------------
  void setDetail(int ures, int vres){
    uStep = 1f/ures;
    vStep = 1f/vres;
  }
  //--------------
  // sets the anchor (corner) points of the patch
  // AND distributes the control points evenly between
  void setCorners(float x0,float y0,float z0,float x1,float y1,float z1,float x2,float y2,float z2,float x3,float y3,float z3){
    patch[0][0] = new float[]{x0,y0,z0};
    patch[1][0] = new float[]{x1,y1,z1};
    patch[0][1] = new float[]{x2,y2,z2};
    patch[1][1] = new float[]{x3,y3,z3};
    interpolateControlPoints();
  }
  //--------------
  // distributes the control points evenly between the anchor points current positions
  void interpolateControlPoints(){
    for(int i=0 ; i<4 ; i++ )
      for(int j=0 ; j<4 ; j++ )
        for(int c=0 ; c<3 ; c++ ){
          float a = ( (3-i)*patch[0][0][c] + i*patch[3][0][c] );
          float b = ( (3-i)*patch[0][3][c] + i*patch[3][3][c] );
          patch[i][j][c] = ( (3-j)*a + j*b ) / 9;
        }
  }
  //--------------
  void setVertex(int i, int j, float x, float y, float z){
    patch[i][j] = new float[]{x,y,z};
  }
  //--------------
  float[][][] getPatch(){
    return patch;
  }
  //--------------
  float x(float u, float v){
    return bezierToCoord( patch, u, v, 0);
  }
  //--------------
  float y(float u, float v){
    return bezierToCoord( patch, u, v, 1);
  }
  //--------------
  float z(float u, float v){
    return bezierToCoord( patch, u, v, 2);
  }
  //--------------
  float normalX(float u, float v){
    float ay = tangentU( patch, u, v, 1 );
    float az = tangentU( patch, u, v, 2 );
    float by = tangentV( patch, u, v, 1 );
    float bz = tangentV( patch, u, v, 2 );
    return ay*bz - az*by;
  }
  //--------------
  float normalY(float u, float v){
    float ax = tangentU( patch, u, v, 0 );
    float az = tangentU( patch, u, v, 2 );
    float bx = tangentV( patch, u, v, 0 );
    float bz = tangentV( patch, u, v, 2 );
    return -ax*bz + az*bx;
  }
  //--------------
  float normalZ(float u, float v){
    float ax = tangentU( patch, u, v, 0 );
    float ay = tangentU( patch, u, v, 1 );
    float bx = tangentV( patch, u, v, 0 );
    float by = tangentV( patch, u, v, 1 );
    return ax*by - ay*bx;
  }
  //--------------
  // return the i:coordinate on a bezierpatch where the patch-coords are (u,v)
  float bezierToCoord(float[][][] patch, float u, float v, int i){
    float u3 = u*u*u;
    float u3p = (1-u)*(1-u)*(1-u);
    float u2 = u*u;
    float u2p = (1-u)*(1-u);
    float v3 = v*v*v;
    float v2 = v*v;
    float v3p = (1-v)*(1-v)*(1-v);
    float v2p = (1-v)*(1-v);

    float a = u3*patch[0][0][i] + 3*u2*(1-u)*patch[1][0][i] + 3*u*u2p*patch[2][0][i] + u3p*patch[3][0][i];
    float b = u3*patch[0][1][i] + 3*u2*(1-u)*patch[1][1][i] + 3*u*u2p*patch[2][1][i] + u3p*patch[3][1][i];
    float c = u3*patch[0][2][i] + 3*u2*(1-u)*patch[1][2][i] + 3*u*u2p*patch[2][2][i] + u3p*patch[3][2][i];
    float d = u3*patch[0][3][i] + 3*u2*(1-u)*patch[1][3][i] + 3*u*u2p*patch[2][3][i] + u3p*patch[3][3][i];

    return v3*a + 3*v2*(1-v)*b + 3*v*v2p*c + v3p*d;
  }
  //--------------
  // return the i:coordinate on a bezierpatch where the patch-coords are (u,v)
  float tangentU(float[][][] patch, float u, float v, int i){
    float u3 = 3*u*u;
    float u3p = -3*(1-u)*(1-u);
    float u2 = 2*u;
    float u2p = -2*(1-u);
    float v3 = v*v*v;
    float v2 = v*v;
    float v3p = (1-v)*(1-v)*(1-v);
    float v2p = (1-v)*(1-v);

    float a = u3*patch[0][0][i] + 3*u2*(1-u)*patch[1][0][i] + 3*u*u2p*patch[2][0][i] + u3p*patch[3][0][i];
    float b = u3*patch[0][1][i] + 3*u2*(1-u)*patch[1][1][i] + 3*u*u2p*patch[2][1][i] + u3p*patch[3][1][i];
    float c = u3*patch[0][2][i] + 3*u2*(1-u)*patch[1][2][i] + 3*u*u2p*patch[2][2][i] + u3p*patch[3][2][i];
    float d = u3*patch[0][3][i] + 3*u2*(1-u)*patch[1][3][i] + 3*u*u2p*patch[2][3][i] + u3p*patch[3][3][i];

    return v3*a + 3*v2*(1-v)*b + 3*v*v2p*c + v3p*d;
  }
  //--------------
  // return the i:coordinate on a bezierpatch where the patch-coords are (u,v)
  float tangentV(float[][][] patch, float u, float v, int i){
    float u3 = u*u*u;
    float u3p = (1-u)*(1-u)*(1-u);
    float u2 = u*u;
    float u2p = (1-u)*(1-u);
    float v3 = 3*v*v;
    float v2 = 2*v;
    float v3p = -3*(1-v)*(1-v);
    float v2p = -2*(1-v);

    float a = u3*patch[0][0][i] + 3*u2*(1-u)*patch[1][0][i] + 3*u*u2p*patch[2][0][i] + u3p*patch[3][0][i];
    float b = u3*patch[0][1][i] + 3*u2*(1-u)*patch[1][1][i] + 3*u*u2p*patch[2][1][i] + u3p*patch[3][1][i];
    float c = u3*patch[0][2][i] + 3*u2*(1-u)*patch[1][2][i] + 3*u*u2p*patch[2][2][i] + u3p*patch[3][2][i];
    float d = u3*patch[0][3][i] + 3*u2*(1-u)*patch[1][3][i] + 3*u*u2p*patch[2][3][i] + u3p*patch[3][3][i];

    return v3*a + 3*v2*(1-v)*b + 3*v*v2p*c + v3p*d;
  }
  //--------------
  // when supplying a PImage as a parameter drawPatch uses the image as a texture
  void drawPatch(PImage skin){
    textureMode(NORMALIZED);

    float u_,v_;

    for(float u=0;u<1;u+=uStep)
      for(float v=0;v<1;v+=vStep){
        u_ = u+uStep;
        v_ = v+vStep;
        beginShape(TRIANGLE_STRIP);
        texture(skin);
        vertex( x(u,v), y(u,v), z(u,v), u,v);
        vertex( x(u_,v), y(u_,v), z(u_,v), u_,v);
        vertex( x(u,v_), y(u,v_), z(u,v_), u,v_);
        vertex( x(u_,v_), y(u_,v_), z(u_,v_), u_,v_);
        endShape(CLOSE);
      }
  }
  //--------------
  // without parameters drawPatch obeys processings stroke and fill settings
  void drawPatch(){
    float u_,v_;
    for(float u=0;u<1;u+=uStep)
      for(float v=0;v<1;v+=vStep){
        u_ = u+uStep;
        v_ = v+vStep;
        beginShape(TRIANGLE_STRIP);
        normal( normalX(u,v), normalY(u,v), normalZ(u,v) );
        vertex( x(u,v), y(u,v), z(u,v));
        normal( normalX(u_,v), normalY(u_,v), normalZ(u_,v) );
        vertex( x(u_,v), y(u_,v), z(u_,v));
        normal( normalX(u,v_), normalY(u,v_), normalZ(u,v_) );
        vertex( x(u,v_), y(u,v_), z(u,v_));
        normal( normalX(u_,v_), normalY(u_,v_), normalZ(u_,v_) );
        vertex( x(u_,v_), y(u_,v_), z(u_,v_));
        endShape(CLOSE);
      }
  }
  //--------------
}