public class cBezier
{
	static int       n = 0; 
	static float[][] bezierMatrix      = null;
	static float[][] bezierJacobian      = null;
	static int[]     splineIndexValues = null;
	static float[]   splineSTValues    = null;
	
	static int[] lods = { 4, 7, 13, 25, 49, 97 };
	
	static cVector bufU = new cVector();
	static cVector bufV = new cVector();
	
    static public BoundaryRep GenPatch( float[] splineControlPoints, int depth, float minu, float maxu, float minv, float maxv )
    {
		int lod = lods[depth];
		
		if ( n != lod )
		{
			n = lod;
			bezierMatrix      = generateBezierMatrix( n );
			bezierJacobian    = generateBezierJacobian( n );
			splineIndexValues = generateSplineIndexValues( n );
			splineSTValues    = generateSplineSTValues( n );
			
			assert( splineIndexValues.length == (n-1)*(n-1)*2*3 );
		}

		//assert( splineControlPoints.length == n*n*3 );

		float[] splinePatchPoints = calculatePatches( splineControlPoints, bezierMatrix, bezierMatrix );  
		float[] splinePatchTanU = calculatePatches( splineControlPoints, bezierJacobian, bezierMatrix );  
		float[] splinePatchTanV = calculatePatches( splineControlPoints, bezierMatrix, bezierJacobian );  

		//assert( splinePatchPoints.length >= n*n*3 );
		
		BoundaryRep brep = new BoundaryRep(n*n, (n-1)*(n-1)*2);
		
		Vertex v;
		
		for( int i=0; i<n*n; i++)
		{
			int i3 = i*3;
		//System.out.println("splinePatchPoints = " + splinePatchPoints[i3] + " " + splinePatchPoints[i3+1] + " " + splinePatchPoints[i3+2]);
			v = new Vertex(splinePatchPoints[i3], splinePatchPoints[i3+1], splinePatchPoints[i3+2]);
			bufU.set(splinePatchTanU[i3], splinePatchTanU[i3+1], splinePatchTanU[i3+2]);
			bufV.set(splinePatchTanV[i3], splinePatchTanV[i3+1], splinePatchTanV[i3+2]);
			v.norm = new cVector();
			v.norm.cross(bufV,bufU);
			v.norm.normalize();
			int i2 = i*2;
			v.s = splineSTValues[i2]*(maxu-minu) + minu;
			v.t = splineSTValues[i2+1]*(maxv-minv) + minv;
			brep.AddVertex(v);
		}
		
		for( int i=0; i<(n-1)*(n-1)*2; i++)
		{
			int i3 = i*3;
		//System.out.println("splineIndexValues = " + splineIndexValues[i3] + " " + splineIndexValues[i3+1] + " " + splineIndexValues[i3+2]);
			brep.AddFace(splineIndexValues[i3], splineIndexValues[i3+1], splineIndexValues[i3+2]);
		}
		
		return brep;
		//setVRMLSplineCoordPoints( splinePatchPoints, splineIndexValues );                        
	}
	
	/***************************  
	*    Bezier Spline math   *  
	***************************/
	
	static float[][] bufferPatches = new float[3][];
	static int patchcount = 0;
	
	static public float[] calculatePatches( float[] p, float[][] bu, float[][] bv) // 0, float[] b1,	float[] b2, float[] b3 )
	{                                               
		int		n  = bu[0].length;   
		
		if (bufferPatches[patchcount] == null || 3*n*n > bufferPatches[patchcount].length)
			bufferPatches[patchcount] = new float[3 * n * n];
		
		float[] pt = bufferPatches[patchcount];
		
		patchcount++;
		patchcount %= 3;
		
		//  println( "calculating Patches" );
		for( int j, k = j = 0; j < n; j++ )
		{
		//System.out.print("b0[" + j + "] = " + b0[j] + ", ");
		//System.out.print("b1[" + j + "] = " + b1[j] + ", ");
		//System.out.print("b2[" + j + "] = " + b2[j] + ", ");
		//System.out.println("b3[" + j + "] = " + b3[j] + ";");
			for( int i = 0;  i < n; i++ )
			{
				pt[k++] = (float)(
					( p[ 0]*bu[0][i] + p[ 3]*bu[1][i] + p[ 6]*bu[2][i] + p[ 9]*bu[3][i]) * bv[0][j] +
					( p[12]*bu[0][i] + p[15]*bu[1][i] + p[18]*bu[2][i] + p[21]*bu[3][i]) * bv[1][j] +
					( p[24]*bu[0][i] + p[27]*bu[1][i] + p[30]*bu[2][i] + p[33]*bu[3][i]) * bv[2][j] +
					( p[36]*bu[0][i] + p[39]*bu[1][i] + p[42]*bu[2][i] + p[45]*bu[3][i]) * bv[3][j] );
				pt[k++] = (float)(
					( p[ 1]*bu[0][i] + p[ 4]*bu[1][i] + p[ 7]*bu[2][i] + p[10]*bu[3][i]) * bv[0][j] +
					( p[13]*bu[0][i] + p[16]*bu[1][i] + p[19]*bu[2][i] + p[22]*bu[3][i]) * bv[1][j] +
					( p[25]*bu[0][i] + p[28]*bu[1][i] + p[31]*bu[2][i] + p[34]*bu[3][i]) * bv[2][j] +
					( p[37]*bu[0][i] + p[40]*bu[1][i] + p[43]*bu[2][i] + p[46]*bu[3][i]) * bv[3][j] );
				pt[k++] = (float)(
					( p[ 2]*bu[0][i] + p[ 5]*bu[1][i] + p[ 8]*bu[2][i] + p[11]*bu[3][i]) * bv[0][j] +
					( p[14]*bu[0][i] + p[17]*bu[1][i] + p[20]*bu[2][i] + p[23]*bu[3][i]) * bv[1][j] +
					( p[26]*bu[0][i] + p[29]*bu[1][i] + p[32]*bu[2][i] + p[35]*bu[3][i]) * bv[2][j] +
					( p[38]*bu[0][i] + p[41]*bu[1][i] + p[44]*bu[2][i] + p[47]*bu[3][i]) * bv[3][j] );
			}                                                                    
		}
		
		return pt;
	}
	
	/*
	static float[] calculatePatches( float[] p, float[][] bezierMatrix )
	{
		return 
		calculatePatches( p, bezierMatrix[0],  bezierMatrix[1],  bezierMatrix[2], bezierMatrix[3] );
	}    
	*/	
	
	
	static public float[][] generateBezierMatrix( int n )
	{
		float          u, u1, u12, u2, u3;                  // original code used all doubles here
		float[][]   valMatrix = new float[4][n];     // 
		float          st = (float) ( 1.0 / ( n - 1 ) );
		
		for( int i = 0; i < n; i++ )
		{
			// Bezier
			u  = i*st; u2 = u*u;
			u1 = 1-u; u12 = u1*u1;
			
			valMatrix[0][i] =     u1 * u12; 
			valMatrix[1][i] = 3 * u  * u12; 
			valMatrix[2][i] = 3 * u2 * u1; 
			//valMatrix[3][i] =     u  * u2;                    

			// Hermite?
			//valMatrix[0][i] = 2*u*u*u - 3*u*u + 1;
			//valMatrix[1][i] = -2*u*u*u + 3*u*u;
			//valMatrix[2][i] = u*u*u - 2*u*u + u;
			//valMatrix[3][i] = u*u*u - u*u;

			/**/
			// B-spline
			u3 = u2*u;
			valMatrix[0][i] = u3/6;
			valMatrix[1][i] = (1 + 3*(u + u2 - u3))/6;
			valMatrix[2][i] = (4 - 6*u2 + 3*u3)/6;
			//valMatrix[3][i] = (1 - 3*(u - u2) - u3)/6;
			/**/

			// random
			//valMatrix[0][i] = u2/5;
			//valMatrix[1][i] = (2 + 1*(u3 + u - u2))/6;
			//valMatrix[2][i] = (4 - 6*u + 4*u3)/3;

			valMatrix[3][i] = 1 - (valMatrix[0][i] + valMatrix[1][i] + valMatrix[2][i]);
		}
		return valMatrix;
	}           //  end generateBezierMatrix method        
	
	static public float[][] generateBezierJacobian( int n )
	{
		float          u, u1, u12, u2;                  // original code used all doubles here
		float[][]   valMatrix = new float[4][n];     // 
		float          st = (float) ( 1.0 / ( n - 1 ) );
		
		for( int i = 0; i < n; i++ )
		{
			u  = i * st; u2  = u*u;
			u1 = 1 - u;  u12 = u1*u1;
			
			// Bezier
			valMatrix[0][i] = -3 * u12; 
			valMatrix[1][i] =  3 * (u12 - u*2*u1); 
			valMatrix[2][i] =  3 * (2*u*u1 - u2); 
			valMatrix[3][i] =  3 * u2;
			
			/**/
			// B_spline
			valMatrix[0][i] = 3 * u2 / 6; 
			valMatrix[1][i] = (3 + 6*u - 9*u2)/6; 
			valMatrix[2][i] = (-12*u + 9*u2)/6; 
			valMatrix[3][i] = (-3 + 6*u - 3*u2)/6;
			/**/
		}
		
		return valMatrix;
	}        
	
	
	static public int[] generateSplineIndexValues( int n )
	{
		int                  NumTriangles       = 2 * (n-1) * (n-1);
		int                  NumVertexValues    = 6 * (n-1) * (n-1);
		int                  NumVertexEntries   = 6 * (n-1) * (n-1);
		
		// Level of Detail : 4 ; 72
		
		//println( "Level of Detail : " + n + " ;  ( " +  NumTriangles        + ", " + NumVertexValues     + ", " + NumVertexEntries    + " ) : " + "  =  2*(n-1)^2 Triangles, 6*(n-1)^2 Values,  8*(n^2) Entries " );
		
		int[] s = new int[NumVertexEntries];
		
		for( int j, k=j=0; j < n-1; j++ )
		{
			for( int i=0;  i < n-1; i++ )
			{
				s[k++] = i+1 + n*j; 
				s[k++] = i   + n*j; 
				s[k++] = i   + n*(j+1); 
				s[k++] = i   + n*(j+1); 
				s[k++] = i+1 + n*(j+1); 
				s[k++] = i+1 + n*j; 
			}
		}   
		return s;
	}

	static public float[] generateSplineSTValues( int n )
	{
		int NumVertexEntries = 2 * n*n;
		
		float step = 1.0f/(n-1);
		
		float[] s = new float[NumVertexEntries];
		
		for( int j=0, k=0; j < n; j++ )
		{
			for( int i=0;  i < n; i++ )
			{
				s[k++] = 1 - i * step; 
				s[k++] = j * step; 
			}
		}   
		return s;
	}
} // class Bezier