From: megacz <adam@megacz.com>
Date: Sat, 28 Feb 2009 21:36:00 +0000 (-0800)
Subject: add first draft of Marching Cubes codefg
X-Git-Url: http://git.megacz.com/?p=anneal.git;a=commitdiff_plain;h=1cff1b9494433864546e120b5fe4952e15b3ea94

add first draft of Marching Cubes codefg
---

diff --git a/src/edu/berkeley/qfat/voxel/MarchingCubes.java b/src/edu/berkeley/qfat/voxel/MarchingCubes.java
new file mode 100644
index 0000000..f8fb4e4
--- /dev/null
+++ b/src/edu/berkeley/qfat/voxel/MarchingCubes.java
@@ -0,0 +1,1024 @@
+package edu.berkeley.qfat.voxel;
+import edu.berkeley.qfat.*;
+import edu.berkeley.qfat.geom.*;
+
+/////////////////////////////////////////////////////////////////////////////////////
+// Derived from:                                                                   //
+// http://local.wasp.uwa.edu.au/~pbourke/geometry/polygonise/marchingsource.cpp    //
+// (which is public domain)                                                        //
+/////////////////////////////////////////////////////////////////////////////////////
+//
+// Marching Cubes Example Program 
+// by Cory Bloyd (corysama@yahoo.com)
+//
+// A simple, portable and complete implementation of the Marching Cubes
+// and Marching Tetrahedrons algorithms in a single source file.
+// There are many ways that this code could be made faster, but the 
+// intent is for the code to be easy to understand.
+//
+// For a description of the algorithm go to
+// http://astronomy.swin.edu.au/pbourke/modelling/polygonise/
+//
+// This code is public domain.
+//
+//////////////////////////////////////////////////////////////////////////////////////
+
+public abstract class MarchingCubes {
+
+    private int iEdgeFlags;
+
+    private static class GLvector {
+        public double fX;
+        public double fY;
+        public double fZ;     
+    };
+
+    // These tables are used so that everything can be done in little
+    // loops that you can look at all at once rather than in pages and
+    // pages of unrolled code.
+    
+    // a2fVertexOffset lists the positions, relative to vertex0, of
+    // each of the 8 vertices of a cube
+    static final double a2fVertexOffset[][] = {
+        {0.0, 0.0, 0.0},{1.0, 0.0, 0.0},{1.0, 1.0, 0.0},{0.0, 1.0, 0.0},
+        {0.0, 0.0, 1.0},{1.0, 0.0, 1.0},{1.0, 1.0, 1.0},{0.0, 1.0, 1.0}
+    };
+    
+    //a2iEdgeConnection lists the index of the endpoint vertices for each of the 12 edges of the cube
+    static final int a2iEdgeConnection[][] = {
+        {0,1}, {1,2}, {2,3}, {3,0},
+        {4,5}, {5,6}, {6,7}, {7,4},
+        {0,4}, {1,5}, {2,6}, {3,7}
+    };
+
+    //a2fEdgeDirection lists the direction vector (vertex1-vertex0) for each edge in the cube
+    static final double a2fEdgeDirection[][] = {
+        {1.0, 0.0, 0.0},{0.0, 1.0, 0.0},{-1.0, 0.0, 0.0},{0.0, -1.0, 0.0},
+        {1.0, 0.0, 0.0},{0.0, 1.0, 0.0},{-1.0, 0.0, 0.0},{0.0, -1.0, 0.0},
+        {0.0, 0.0, 1.0},{0.0, 0.0, 1.0},{ 0.0, 0.0, 1.0},{0.0,  0.0, 1.0}
+    };
+
+    // a2iTetrahedronEdgeConnection lists the index of the endpoint
+    // vertices for each of the 6 edges of the tetrahedron
+    static final int a2iTetrahedronEdgeConnection[][] = {
+        {0,1},  {1,2},  {2,0},  {0,3},  {1,3},  {2,3}
+    };
+
+    // a2iTetrahedronEdgeConnection lists the index of verticies from a cube 
+    // that made up each of the six tetrahedrons within the cube
+    static final int a2iTetrahedronsInACube[][] = {
+        {0,5,1,6},
+        {0,1,2,6},
+        {0,2,3,6},
+        {0,3,7,6},
+        {0,7,4,6},
+        {0,4,5,6},
+    };
+
+    static final double afAmbientWhite [] = {0.25, 0.25, 0.25, 1.00}; 
+    static final double afAmbientRed   [] = {0.25, 0.00, 0.00, 1.00}; 
+    static final double afAmbientGreen [] = {0.00, 0.25, 0.00, 1.00}; 
+    static final double afAmbientBlue  [] = {0.00, 0.00, 0.25, 1.00}; 
+    static final double afDiffuseWhite [] = {0.75, 0.75, 0.75, 1.00}; 
+    static final double afDiffuseRed   [] = {0.75, 0.00, 0.00, 1.00}; 
+    static final double afDiffuseGreen [] = {0.00, 0.75, 0.00, 1.00}; 
+    static final double afDiffuseBlue  [] = {0.00, 0.00, 0.75, 1.00}; 
+    static final double afSpecularWhite[] = {1.00, 1.00, 1.00, 1.00}; 
+    static final double afSpecularRed  [] = {1.00, 0.25, 0.25, 1.00}; 
+    static final double afSpecularGreen[] = {0.25, 1.00, 0.25, 1.00}; 
+    static final double afSpecularBlue [] = {0.25, 0.25, 1.00, 1.00}; 
+
+    GLvector  sSourcePoint[] = new GLvector[3];
+
+    /*
+      GLenum    ePolygonMode = GL_FILL;
+      int     iDataSetSize = 16;
+      double   fStepSize = 1.0/iDataSetSize;
+      double   fTargetValue = 48.0;
+      double   fTime = 0.0;
+      boolean bSpin = true;
+      boolean bMove = true;
+      boolean bLight = true;
+
+
+      void vIdle();
+      void vDrawScene(); 
+      void vResize(GLsizei, GLsizei);
+      void vKeyboard(unsigned char cKey, int iX, int iY);
+      void vSpecial(int iKey, int iX, int iY);
+
+      GLvoid vPrintHelp();
+      GLvoid vSetTime(double fTime);
+      double fSample1(double fX, double fY, double fZ);
+      double fSample2(double fX, double fY, double fZ);
+      double fSample3(double fX, double fY, double fZ);
+      double (*fSample)(double fX, double fY, double fZ) = fSample1;
+
+      GLvoid vMarchingCubes();
+      GLvoid vMarchCube1(double fX, double fY, double fZ, double fScale);
+      GLvoid vMarchCube2(double fX, double fY, double fZ, double fScale);
+      GLvoid (*vMarchCube)(double fX, double fY, double fZ, double fScale) = vMarchCube1;
+
+      void main(int argc, char **argv) 
+      { 
+      double afPropertiesAmbient [] = {0.50, 0.50, 0.50, 1.00}; 
+      double afPropertiesDiffuse [] = {0.75, 0.75, 0.75, 1.00}; 
+      double afPropertiesSpecular[] = {1.00, 1.00, 1.00, 1.00}; 
+
+      GLsizei iWidth = 640.0; 
+      GLsizei iHeight = 480.0; 
+
+      glutInit(&argc, argv);
+      glutInitWindowPosition( 0, 0);
+      glutInitWindowSize(iWidth, iHeight);
+      glutInitDisplayMode( GLUT_RGB | GLUT_DEPTH | GLUT_DOUBLE );
+      glutCreateWindow( "Marching Cubes" );
+      glutDisplayFunc( vDrawScene );
+      glutIdleFunc( vIdle );
+      glutReshapeFunc( vResize );
+      glutKeyboardFunc( vKeyboard );
+      glutSpecialFunc( vSpecial );
+
+      glClearColor( 0.0, 0.0, 0.0, 1.0 ); 
+      glClearDepth( 1.0 ); 
+
+      glEnable(GL_DEPTH_TEST); 
+      glEnable(GL_LIGHTING);
+      glPolygonMode(GL_FRONT_AND_BACK, ePolygonMode);
+
+      glLightfv( GL_LIGHT0, GL_AMBIENT,  afPropertiesAmbient); 
+      glLightfv( GL_LIGHT0, GL_DIFFUSE,  afPropertiesDiffuse); 
+      glLightfv( GL_LIGHT0, GL_SPECULAR, afPropertiesSpecular); 
+      glLightModelf(GL_LIGHT_MODEL_TWO_SIDE, 1.0); 
+
+      glEnable( GL_LIGHT0 ); 
+
+      glMaterialfv(GL_BACK,  GL_AMBIENT,   afAmbientGreen); 
+      glMaterialfv(GL_BACK,  GL_DIFFUSE,   afDiffuseGreen); 
+      glMaterialfv(GL_FRONT, GL_AMBIENT,   afAmbientBlue); 
+      glMaterialfv(GL_FRONT, GL_DIFFUSE,   afDiffuseBlue); 
+      glMaterialfv(GL_FRONT, GL_SPECULAR,  afSpecularWhite); 
+      glMaterialf( GL_FRONT, GL_SHININESS, 25.0); 
+
+      vResize(iWidth, iHeight); 
+
+      vPrintHelp();
+      glutMainLoop(); 
+      }
+
+      GLvoid vPrintHelp()
+      {
+      printf("Marching Cubes Example by Cory Bloyd (dejaspaminacan@my-deja.com)\n\n");
+
+      printf("+/-  increase/decrease sample density\n");
+      printf("PageUp/PageDown  increase/decrease surface value\n");
+      printf("s  change sample function\n");
+      printf("c  toggle marching cubes / marching tetrahedrons\n");
+      printf("w  wireframe on/off\n");
+      printf("l  toggle lighting / color-by-normal\n");
+      printf("Home  spin scene on/off\n");
+      printf("End  source point animation on/off\n");
+      }
+
+
+      void vResize( GLsizei iWidth, GLsizei iHeight ) 
+      { 
+      double fAspect, fHalfWorldSize = (1.4142135623730950488016887242097/2); 
+
+      glViewport( 0, 0, iWidth, iHeight ); 
+      glMatrixMode (GL_PROJECTION);
+      glLoadIdentity ();
+
+      if(iWidth <= iHeight)
+      {
+      fAspect = (double)iHeight / (double)iWidth;
+      glOrtho(-fHalfWorldSize, fHalfWorldSize, -fHalfWorldSize*fAspect,
+      fHalfWorldSize*fAspect, -10*fHalfWorldSize, 10*fHalfWorldSize);
+      }
+      else
+      {
+      fAspect = (double)iWidth / (double)iHeight; 
+      glOrtho(-fHalfWorldSize*fAspect, fHalfWorldSize*fAspect, -fHalfWorldSize,
+      fHalfWorldSize, -10*fHalfWorldSize, 10*fHalfWorldSize);
+      }
+ 
+      glMatrixMode( GL_MODELVIEW );
+      }
+
+      void vKeyboard(unsigned char cKey, int iX, int iY)
+      {
+      switch(cKey)
+      {
+      case 'w' :
+      {
+      if(ePolygonMode == GL_LINE)
+      {
+      ePolygonMode = GL_FILL;
+      }
+      else
+      {
+      ePolygonMode = GL_LINE;
+      }
+      glPolygonMode(GL_FRONT_AND_BACK, ePolygonMode);
+      } break;
+      case '+' :
+      case '=' :
+      {
+      ++iDataSetSize;
+      fStepSize = 1.0/iDataSetSize;
+      } break;
+      case '-' :
+      {
+      if(iDataSetSize > 1)
+      {
+      --iDataSetSize;
+      fStepSize = 1.0/iDataSetSize;
+      }
+      } break;
+      case 'c' :
+      {
+      if(vMarchCube == vMarchCube1)
+      {
+      vMarchCube = vMarchCube2;//Use Marching Tetrahedrons
+      }
+      else
+      {
+      vMarchCube = vMarchCube1;//Use Marching Cubes
+      }
+      } break;
+      case 's' :
+      {
+      if(fSample == fSample1)
+      {
+      fSample = fSample2;
+      }
+      else if(fSample == fSample2)
+      {
+      fSample = fSample3;
+      }
+      else
+      {
+      fSample = fSample1;
+      }
+      } break;
+      case 'l' :
+      {
+      if(bLight)
+      {
+      glDisable(GL_LIGHTING);//use vertex colors
+      }
+      else
+      {
+      glEnable(GL_LIGHTING);//use lit material color
+      }
+
+      bLight = !bLight;
+      };
+      }
+      }
+
+
+      void vSpecial(int iKey, int iX, int iY)
+      {
+      switch(iKey)
+      {
+      case GLUT_KEY_PAGE_UP :
+      {
+      if(fTargetValue < 1000.0)
+      {
+      fTargetValue *= 1.1;
+      }
+      } break;
+      case GLUT_KEY_PAGE_DOWN :
+      {
+      if(fTargetValue > 1.0)
+      {
+      fTargetValue /= 1.1;
+      }
+      } break;
+      case GLUT_KEY_HOME :
+      {
+      bSpin = !bSpin;
+      } break;
+      case GLUT_KEY_END :
+      {
+      bMove = !bMove;
+      } break;
+      }
+      }
+
+      void vIdle()
+      {
+      glutPostRedisplay();
+      }
+
+      void vDrawScene() 
+      { 
+      static double fPitch = 0.0;
+      static double fYaw   = 0.0;
+      static double fTime = 0.0;
+
+      glClear( GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT ); 
+
+      glPushMatrix(); 
+
+      if(bSpin)
+      {
+      fPitch += 4.0;
+      fYaw   += 2.5;
+      }
+      if(bMove)
+      {
+      fTime  += 0.025;
+      }
+
+      vSetTime(fTime);
+
+      glTranslatef(0.0, 0.0, -1.0);  
+      glRotatef( -fPitch, 1.0, 0.0, 0.0);
+      glRotatef(     0.0, 0.0, 1.0, 0.0);
+      glRotatef(    fYaw, 0.0, 0.0, 1.0);
+
+      glPushAttrib(GL_LIGHTING_BIT);
+      glDisable(GL_LIGHTING);
+      glColor3f(1.0, 1.0, 1.0);
+      glutWireCube(1.0); 
+      glPopAttrib(); 
+
+
+      glPushMatrix(); 
+      glTranslatef(-0.5, -0.5, -0.5);
+      glBegin(GL_TRIANGLES);
+      vMarchingCubes();
+      glEnd();
+      glPopMatrix(); 
+
+
+      glPopMatrix(); 
+
+      glutSwapBuffers(); 
+      }
+
+
+
+      //vGetColor generates a color from a given position and normal of a point
+      GLvoid vGetColor(GLvector &rfColor, GLvector &rfPosition, GLvector &rfNormal)
+      {
+      double fX = rfNormal.fX;
+      double fY = rfNormal.fY;
+      double fZ = rfNormal.fZ;
+      rfColor.fX = (fX > 0.0 ? fX : 0.0) + (fY < 0.0 ? -0.5*fY : 0.0) + (fZ < 0.0 ? -0.5*fZ : 0.0);
+      rfColor.fY = (fY > 0.0 ? fY : 0.0) + (fZ < 0.0 ? -0.5*fZ : 0.0) + (fX < 0.0 ? -0.5*fX : 0.0);
+      rfColor.fZ = (fZ > 0.0 ? fZ : 0.0) + (fX < 0.0 ? -0.5*fX : 0.0) + (fY < 0.0 ? -0.5*fY : 0.0);
+      }
+
+
+
+      //Generate a sample data set.  fSample1(), fSample2() and fSample3() define three scalar fields whose
+      // values vary by the X,Y and Z coordinates and by the fTime value set by vSetTime()
+      GLvoid vSetTime(double fNewTime)
+      {
+      double fOffset;
+      int iSourceNum;
+
+      for(iSourceNum = 0; iSourceNum < 3; iSourceNum++)
+      {
+      sSourcePoint[iSourceNum].fX = 0.5;
+      sSourcePoint[iSourceNum].fY = 0.5;
+      sSourcePoint[iSourceNum].fZ = 0.5;
+      }
+
+      fTime = fNewTime;
+      fOffset = 1.0 + sinf(fTime);
+      sSourcePoint[0].fX *= fOffset;
+      sSourcePoint[1].fY *= fOffset;
+      sSourcePoint[2].fZ *= fOffset;
+      }
+
+
+      //fSample1 finds the distance of (fX, fY, fZ) from three moving points
+      double fSample1(double fX, double fY, double fZ)
+      {
+      GLdouble fResult = 0.0;
+      GLdouble fDx, fDy, fDz;
+      fDx = fX - sSourcePoint[0].fX;
+      fDy = fY - sSourcePoint[0].fY;
+      fDz = fZ - sSourcePoint[0].fZ;
+      fResult += 0.5/(fDx*fDx + fDy*fDy + fDz*fDz);
+
+      fDx = fX - sSourcePoint[1].fX;
+      fDy = fY - sSourcePoint[1].fY;
+      fDz = fZ - sSourcePoint[1].fZ;
+      fResult += 1.0/(fDx*fDx + fDy*fDy + fDz*fDz);
+
+      fDx = fX - sSourcePoint[2].fX;
+      fDy = fY - sSourcePoint[2].fY;
+      fDz = fZ - sSourcePoint[2].fZ;
+      fResult += 1.5/(fDx*fDx + fDy*fDy + fDz*fDz);
+
+      return fResult;
+      }
+
+      //fSample2 finds the distance of (fX, fY, fZ) from three moving lines
+      double fSample2(double fX, double fY, double fZ)
+      {
+      GLdouble fResult = 0.0;
+      GLdouble fDx, fDy, fDz;
+      fDx = fX - sSourcePoint[0].fX;
+      fDy = fY - sSourcePoint[0].fY;
+      fResult += 0.5/(fDx*fDx + fDy*fDy);
+
+      fDx = fX - sSourcePoint[1].fX;
+      fDz = fZ - sSourcePoint[1].fZ;
+      fResult += 0.75/(fDx*fDx + fDz*fDz);
+
+      fDy = fY - sSourcePoint[2].fY;
+      fDz = fZ - sSourcePoint[2].fZ;
+      fResult += 1.0/(fDy*fDy + fDz*fDz);
+
+      return fResult;
+      }
+
+
+      //fSample2 defines a height field by plugging the distance from the center into the sin and cos functions
+      double fSample3(double fX, double fY, double fZ)
+      {
+      double fHeight = 20.0*(fTime + sqrt((0.5-fX)*(0.5-fX) + (0.5-fY)*(0.5-fY)));
+      fHeight = 1.5 + 0.1*(sinf(fHeight) + cosf(fHeight));
+      GLdouble fResult = (fHeight - fZ)*50.0;
+
+      return fResult;
+      }
+
+
+
+
+
+
+      //vMarchTetrahedron performs the Marching Tetrahedrons algorithm on a single tetrahedron
+      GLvoid vMarchTetrahedron(GLvector *pasTetrahedronPosition, double *pafTetrahedronValue)
+      {
+      extern int aiTetrahedronEdgeFlags[16];
+      extern int a2iTetrahedronTriangles[16][7];
+
+      int iEdge, iVert0, iVert1, iEdgeFlags, iTriangle, iCorner, iVertex, iFlagIndex = 0;
+      double fOffset, fInvOffset, fValue = 0.0;
+      GLvector asEdgeVertex[6];
+      GLvector asEdgeNorm[6];
+      GLvector sColor;
+
+      //Find which vertices are inside of the surface and which are outside
+      for(iVertex = 0; iVertex < 4; iVertex++)
+      {
+      if(pafTetrahedronValue[iVertex] <= fTargetValue) 
+      iFlagIndex |= 1<<iVertex;
+      }
+
+      //Find which edges are intersected by the surface
+      iEdgeFlags = aiTetrahedronEdgeFlags[iFlagIndex];
+
+      //If the tetrahedron is entirely inside or outside of the surface, then there will be no intersections
+      if(iEdgeFlags == 0)
+      {
+      return;
+      }
+      //Find the point of intersection of the surface with each edge
+      // Then find the normal to the surface at those points
+      for(iEdge = 0; iEdge < 6; iEdge++)
+      {
+      //if there is an intersection on this edge
+      if(iEdgeFlags & (1<<iEdge))
+      {
+      iVert0 = a2iTetrahedronEdgeConnection[iEdge][0];
+      iVert1 = a2iTetrahedronEdgeConnection[iEdge][1];
+      fOffset = fGetOffset(pafTetrahedronValue[iVert0], pafTetrahedronValue[iVert1], fTargetValue);
+      fInvOffset = 1.0 - fOffset;
+
+      asEdgeVertex[iEdge].fX = fInvOffset*pasTetrahedronPosition[iVert0].fX  +  fOffset*pasTetrahedronPosition[iVert1].fX;
+      asEdgeVertex[iEdge].fY = fInvOffset*pasTetrahedronPosition[iVert0].fY  +  fOffset*pasTetrahedronPosition[iVert1].fY;
+      asEdgeVertex[iEdge].fZ = fInvOffset*pasTetrahedronPosition[iVert0].fZ  +  fOffset*pasTetrahedronPosition[iVert1].fZ;
+                        
+      vGetNormal(asEdgeNorm[iEdge], asEdgeVertex[iEdge].fX, asEdgeVertex[iEdge].fY, asEdgeVertex[iEdge].fZ);
+      }
+      }
+      //Draw the triangles that were found.  There can be up to 2 per tetrahedron
+      for(iTriangle = 0; iTriangle < 2; iTriangle++)
+      {
+      if(a2iTetrahedronTriangles[iFlagIndex][3*iTriangle] < 0)
+      break;
+
+      for(iCorner = 0; iCorner < 3; iCorner++)
+      {
+      iVertex = a2iTetrahedronTriangles[iFlagIndex][3*iTriangle+iCorner];
+
+      vGetColor(sColor, asEdgeVertex[iVertex], asEdgeNorm[iVertex]);
+      glColor3f(sColor.fX, sColor.fY, sColor.fZ);
+      glNormal3f(asEdgeNorm[iVertex].fX,   asEdgeNorm[iVertex].fY,   asEdgeNorm[iVertex].fZ);
+      glVertex3f(asEdgeVertex[iVertex].fX, asEdgeVertex[iVertex].fY, asEdgeVertex[iVertex].fZ);
+      }
+      }
+      }
+
+
+
+      //vMarchCube2 performs the Marching Tetrahedrons algorithm on a single cube by making six calls to vMarchTetrahedron
+      GLvoid vMarchCube2(double fX, double fY, double fZ, double fScale)
+      {
+      int iVertex, iTetrahedron, iVertexInACube;
+      GLvector asCubePosition[8];
+      double  afCubeValue[8];
+      GLvector asTetrahedronPosition[4];
+      double  afTetrahedronValue[4];
+
+      //Make a local copy of the cube's corner positions
+      for(iVertex = 0; iVertex < 8; iVertex++)
+      {
+      asCubePosition[iVertex].fX = fX + a2fVertexOffset[iVertex][0]*fScale;
+      asCubePosition[iVertex].fY = fY + a2fVertexOffset[iVertex][1]*fScale;
+      asCubePosition[iVertex].fZ = fZ + a2fVertexOffset[iVertex][2]*fScale;
+      }
+
+      //Make a local copy of the cube's corner values
+      for(iVertex = 0; iVertex < 8; iVertex++)
+      {
+      afCubeValue[iVertex] = fSample(asCubePosition[iVertex].fX,
+      asCubePosition[iVertex].fY,
+      asCubePosition[iVertex].fZ);
+      }
+
+      for(iTetrahedron = 0; iTetrahedron < 6; iTetrahedron++)
+      {
+      for(iVertex = 0; iVertex < 4; iVertex++)
+      {
+      iVertexInACube = a2iTetrahedronsInACube[iTetrahedron][iVertex];
+      asTetrahedronPosition[iVertex].fX = asCubePosition[iVertexInACube].fX;
+      asTetrahedronPosition[iVertex].fY = asCubePosition[iVertexInACube].fY;
+      asTetrahedronPosition[iVertex].fZ = asCubePosition[iVertexInACube].fZ;
+      afTetrahedronValue[iVertex] = afCubeValue[iVertexInACube];
+      }
+      vMarchTetrahedron(asTetrahedronPosition, afTetrahedronValue);
+      }
+      }
+    */        
+
+
+    // vMarchingCubes iterates over the entire dataset, calling vMarchCube on each cube
+    void vMarchingCubes(int iDataSetSize, int fStepSize)
+    {
+        int iX, iY, iZ;
+        for(iX = 0; iX < iDataSetSize; iX++)
+            for(iY = 0; iY < iDataSetSize; iY++)
+                for(iZ = 0; iZ < iDataSetSize; iZ++)
+                    {
+                        vMarchCube(iX*fStepSize, iY*fStepSize, iZ*fStepSize, fStepSize);
+                    }
+    }
+
+    public abstract double fSample(double x, double y, double z);
+    public abstract double fTargetValue();
+
+    //vGetNormal() finds the gradient of the scalar field at a point
+    //This gradient can be used as a very accurate vertx normal for lighting calculations
+    void vGetNormal(GLvector rfNormal, double fX, double fY, double fZ) {
+        rfNormal.fX = fSample(fX-0.01, fY, fZ) - fSample(fX+0.01, fY, fZ);
+        rfNormal.fY = fSample(fX, fY-0.01, fZ) - fSample(fX, fY+0.01, fZ);
+        rfNormal.fZ = fSample(fX, fY, fZ-0.01) - fSample(fX, fY, fZ+0.01);
+        vNormalizeVector(rfNormal, rfNormal);
+    }
+
+    void vNormalizeVector(GLvector rfVectorResult, GLvector rfVectorSource) {
+        double fOldLength;
+        double fScale;
+        
+        fOldLength = Math.sqrt( (rfVectorSource.fX * rfVectorSource.fX) +
+                                (rfVectorSource.fY * rfVectorSource.fY) +
+                                (rfVectorSource.fZ * rfVectorSource.fZ) );
+        
+        if(fOldLength == 0.0) {
+            rfVectorResult.fX = rfVectorSource.fX;
+            rfVectorResult.fY = rfVectorSource.fY;
+            rfVectorResult.fZ = rfVectorSource.fZ;
+        } else {
+            fScale = 1.0/fOldLength;
+            rfVectorResult.fX = rfVectorSource.fX*fScale;
+            rfVectorResult.fY = rfVectorSource.fY*fScale;
+            rfVectorResult.fZ = rfVectorSource.fZ*fScale;
+        }
+    }
+
+    // fGetOffset finds the approximate point of intersection of the surface
+    // between two points with the values fValue1 and fValue2
+    double fGetOffset(double fValue1, double fValue2, double fValueDesired) {
+        double fDelta = fValue2 - fValue1;
+        if(fDelta == 0.0) return 0.5;
+        return (fValueDesired - fValue1)/fDelta;
+    }
+
+    //vMarchCube performs the Marching Cubes algorithm on a single cube
+    void vMarchCube(double fX, double fY, double fZ, double fScale) {
+        int iCorner, iVertex, iVertexTest, iEdge, iTriangle, iFlagIndex, iEdgeFlags;
+        double fOffset;
+        GLvector sColor;
+        double afCubeValue[] = new double[8];
+        GLvector asEdgeVertex[] = new GLvector[12];
+        GLvector asEdgeNorm[] = new GLvector[12];
+
+        //Make a local copy of the values at the cube's corners
+        for(iVertex = 0; iVertex < 8; iVertex++)
+            afCubeValue[iVertex] = fSample(fX + a2fVertexOffset[iVertex][0]*fScale,
+                                           fY + a2fVertexOffset[iVertex][1]*fScale,
+                                           fZ + a2fVertexOffset[iVertex][2]*fScale);
+        
+        //Find which vertices are inside of the surface and which are outside
+        iFlagIndex = 0;
+        for(iVertexTest = 0; iVertexTest < 8; iVertexTest++) {
+            if(afCubeValue[iVertexTest] <= fTargetValue()) 
+                iFlagIndex |= 1<<iVertexTest;
+        }
+        
+        //Find which edges are intersected by the surface
+        iEdgeFlags = aiCubeEdgeFlags[iFlagIndex];
+        
+        //If the cube is entirely inside or outside of the surface, then there will be no intersections
+        if(iEdgeFlags == 0) return;
+
+        //Find the point of intersection of the surface with each edge
+        //Then find the normal to the surface at those points
+        for(iEdge = 0; iEdge < 12; iEdge++) {
+            //if there is an intersection on this edge
+            if ((iEdgeFlags & (1<<iEdge))!=0) {
+                fOffset = fGetOffset(afCubeValue[ a2iEdgeConnection[iEdge][0] ], 
+                                     afCubeValue[ a2iEdgeConnection[iEdge][1] ], fTargetValue());
+                
+                asEdgeVertex[iEdge].fX = fX + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][0]  +  fOffset * a2fEdgeDirection[iEdge][0]) * fScale;
+                asEdgeVertex[iEdge].fY = fY + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][1]  +  fOffset * a2fEdgeDirection[iEdge][1]) * fScale;
+                asEdgeVertex[iEdge].fZ = fZ + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][2]  +  fOffset * a2fEdgeDirection[iEdge][2]) * fScale;
+                
+                vGetNormal(asEdgeNorm[iEdge], asEdgeVertex[iEdge].fX, asEdgeVertex[iEdge].fY, asEdgeVertex[iEdge].fZ);
+            }
+        }
+
+        // Draw the triangles that were found.  There can be up to five per cube
+        for(iTriangle = 0; iTriangle < 5; iTriangle++) {
+            if(a2iTriangleConnectionTable[iFlagIndex][3*iTriangle] < 0)
+                break;
+
+            for(iCorner = 0; iCorner < 3; iCorner++) {
+                iVertex = a2iTriangleConnectionTable[iFlagIndex][3*iTriangle+iCorner];
+                /*
+                vGetColor(sColor, asEdgeVertex[iVertex], asEdgeNorm[iVertex]);
+                glColor3f(sColor.fX, sColor.fY, sColor.fZ);
+                glNormal3f(asEdgeNorm[iVertex].fX,   asEdgeNorm[iVertex].fY,   asEdgeNorm[iVertex].fZ);
+                glVertex3f(asEdgeVertex[iVertex].fX, asEdgeVertex[iVertex].fY, asEdgeVertex[iVertex].fZ);
+                */
+                //mesh.newT(
+            }
+        }
+    }
+
+
+    // For any edge, if one vertex is inside of the surface and the
+    //  other is outside of the surface then the edge intersects the
+    //  surface For each of the 4 vertices of the tetrahedron can be
+    //  two possible states : either inside or outside of the surface
+    //  For any tetrahedron the are 2^4=16 possible sets of vertex
+    //  states This table lists the edges intersected by the surface
+    //  for all 16 possible vertex states There are 6 edges.  For each
+    //  entry in the table, if edge #n is intersected, then bit #n is
+    //  set to 1
+
+    static int aiTetrahedronEdgeFlags[] = {
+        0x00, 0x0d, 0x13, 0x1e, 0x26, 0x2b, 0x35, 0x38, 0x38, 0x35, 0x2b, 0x26, 0x1e, 0x13, 0x0d, 0x00, 
+    };
+
+
+    // For each of the possible vertex states listed in
+    // aiTetrahedronEdgeFlags there is a specific triangulation of the
+    // edge intersection points.  a2iTetrahedronTriangles lists all of
+    // them in the form of 0-2 edge triples with the list terminated
+    // by the invalid value -1.
+    //
+    // I generated this table by hand
+
+    static int a2iTetrahedronTriangles[][] = {
+        {-1, -1, -1, -1, -1, -1, -1},
+        { 0,  3,  2, -1, -1, -1, -1},
+        { 0,  1,  4, -1, -1, -1, -1},
+        { 1,  4,  2,  2,  4,  3, -1},
+
+        { 1,  2,  5, -1, -1, -1, -1},
+        { 0,  3,  5,  0,  5,  1, -1},
+        { 0,  2,  5,  0,  5,  4, -1},
+        { 5,  4,  3, -1, -1, -1, -1},
+
+        { 3,  4,  5, -1, -1, -1, -1},
+        { 4,  5,  0,  5,  2,  0, -1},
+        { 1,  5,  0,  5,  3,  0, -1},
+        { 5,  2,  1, -1, -1, -1, -1},
+
+        { 3,  4,  2,  2,  4,  1, -1},
+        { 4,  1,  0, -1, -1, -1, -1},
+        { 2,  3,  0, -1, -1, -1, -1},
+        {-1, -1, -1, -1, -1, -1, -1},
+    };
+
+    // For any edge, if one vertex is inside of the surface and the
+    //  other is outside of the surface then the edge intersects the
+    //  surface For each of the 8 vertices of the cube can be two
+    //  possible states : either inside or outside of the surface For
+    //  any cube the are 2^8=256 possible sets of vertex states This
+    //  table lists the edges intersected by the surface for all 256
+    //  possible vertex states There are 12 edges.  For each entry in
+    //  the table, if edge #n is intersected, then bit #n is set to 1
+
+    int aiCubeEdgeFlags[] = {
+        0x000, 0x109, 0x203, 0x30a, 0x406, 0x50f, 0x605, 0x70c, 0x80c, 0x905, 0xa0f, 0xb06, 0xc0a, 0xd03, 0xe09, 0xf00, 
+        0x190, 0x099, 0x393, 0x29a, 0x596, 0x49f, 0x795, 0x69c, 0x99c, 0x895, 0xb9f, 0xa96, 0xd9a, 0xc93, 0xf99, 0xe90, 
+        0x230, 0x339, 0x033, 0x13a, 0x636, 0x73f, 0x435, 0x53c, 0xa3c, 0xb35, 0x83f, 0x936, 0xe3a, 0xf33, 0xc39, 0xd30, 
+        0x3a0, 0x2a9, 0x1a3, 0x0aa, 0x7a6, 0x6af, 0x5a5, 0x4ac, 0xbac, 0xaa5, 0x9af, 0x8a6, 0xfaa, 0xea3, 0xda9, 0xca0, 
+        0x460, 0x569, 0x663, 0x76a, 0x066, 0x16f, 0x265, 0x36c, 0xc6c, 0xd65, 0xe6f, 0xf66, 0x86a, 0x963, 0xa69, 0xb60, 
+        0x5f0, 0x4f9, 0x7f3, 0x6fa, 0x1f6, 0x0ff, 0x3f5, 0x2fc, 0xdfc, 0xcf5, 0xfff, 0xef6, 0x9fa, 0x8f3, 0xbf9, 0xaf0, 
+        0x650, 0x759, 0x453, 0x55a, 0x256, 0x35f, 0x055, 0x15c, 0xe5c, 0xf55, 0xc5f, 0xd56, 0xa5a, 0xb53, 0x859, 0x950, 
+        0x7c0, 0x6c9, 0x5c3, 0x4ca, 0x3c6, 0x2cf, 0x1c5, 0x0cc, 0xfcc, 0xec5, 0xdcf, 0xcc6, 0xbca, 0xac3, 0x9c9, 0x8c0, 
+        0x8c0, 0x9c9, 0xac3, 0xbca, 0xcc6, 0xdcf, 0xec5, 0xfcc, 0x0cc, 0x1c5, 0x2cf, 0x3c6, 0x4ca, 0x5c3, 0x6c9, 0x7c0, 
+        0x950, 0x859, 0xb53, 0xa5a, 0xd56, 0xc5f, 0xf55, 0xe5c, 0x15c, 0x055, 0x35f, 0x256, 0x55a, 0x453, 0x759, 0x650, 
+        0xaf0, 0xbf9, 0x8f3, 0x9fa, 0xef6, 0xfff, 0xcf5, 0xdfc, 0x2fc, 0x3f5, 0x0ff, 0x1f6, 0x6fa, 0x7f3, 0x4f9, 0x5f0, 
+        0xb60, 0xa69, 0x963, 0x86a, 0xf66, 0xe6f, 0xd65, 0xc6c, 0x36c, 0x265, 0x16f, 0x066, 0x76a, 0x663, 0x569, 0x460, 
+        0xca0, 0xda9, 0xea3, 0xfaa, 0x8a6, 0x9af, 0xaa5, 0xbac, 0x4ac, 0x5a5, 0x6af, 0x7a6, 0x0aa, 0x1a3, 0x2a9, 0x3a0, 
+        0xd30, 0xc39, 0xf33, 0xe3a, 0x936, 0x83f, 0xb35, 0xa3c, 0x53c, 0x435, 0x73f, 0x636, 0x13a, 0x033, 0x339, 0x230, 
+        0xe90, 0xf99, 0xc93, 0xd9a, 0xa96, 0xb9f, 0x895, 0x99c, 0x69c, 0x795, 0x49f, 0x596, 0x29a, 0x393, 0x099, 0x190, 
+        0xf00, 0xe09, 0xd03, 0xc0a, 0xb06, 0xa0f, 0x905, 0x80c, 0x70c, 0x605, 0x50f, 0x406, 0x30a, 0x203, 0x109, 0x000
+    };
+
+    //  For each of the possible vertex states listed in
+    //  aiCubeEdgeFlags there is a specific triangulation of the edge
+    //  intersection points.  a2iTriangleConnectionTable lists all of
+    //  them in the form of 0-5 edge triples with the list terminated
+    //  by the invalid value -1.  For example:
+    //  a2iTriangleConnectionTable[3] list the 2 triangles formed when
+    //  corner[0] and corner[1] are inside of the surface, but the
+    //  rest of the cube is not.
+    //
+    //  I found this table in an example program someone wrote long
+    //  ago.  It was probably generated by hand
+
+    int a2iTriangleConnectionTable[][] = {
+        {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 8, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 1, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 8, 3, 9, 8, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 8, 3, 1, 2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {9, 2, 10, 0, 2, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {2, 8, 3, 2, 10, 8, 10, 9, 8, -1, -1, -1, -1, -1, -1, -1},
+        {3, 11, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 11, 2, 8, 11, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 9, 0, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 11, 2, 1, 9, 11, 9, 8, 11, -1, -1, -1, -1, -1, -1, -1},
+        {3, 10, 1, 11, 10, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 10, 1, 0, 8, 10, 8, 11, 10, -1, -1, -1, -1, -1, -1, -1},
+        {3, 9, 0, 3, 11, 9, 11, 10, 9, -1, -1, -1, -1, -1, -1, -1},
+        {9, 8, 10, 10, 8, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {4, 7, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {4, 3, 0, 7, 3, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 1, 9, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {4, 1, 9, 4, 7, 1, 7, 3, 1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 2, 10, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {3, 4, 7, 3, 0, 4, 1, 2, 10, -1, -1, -1, -1, -1, -1, -1},
+        {9, 2, 10, 9, 0, 2, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1},
+        {2, 10, 9, 2, 9, 7, 2, 7, 3, 7, 9, 4, -1, -1, -1, -1},
+        {8, 4, 7, 3, 11, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {11, 4, 7, 11, 2, 4, 2, 0, 4, -1, -1, -1, -1, -1, -1, -1},
+        {9, 0, 1, 8, 4, 7, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1},
+        {4, 7, 11, 9, 4, 11, 9, 11, 2, 9, 2, 1, -1, -1, -1, -1},
+        {3, 10, 1, 3, 11, 10, 7, 8, 4, -1, -1, -1, -1, -1, -1, -1},
+        {1, 11, 10, 1, 4, 11, 1, 0, 4, 7, 11, 4, -1, -1, -1, -1},
+        {4, 7, 8, 9, 0, 11, 9, 11, 10, 11, 0, 3, -1, -1, -1, -1},
+        {4, 7, 11, 4, 11, 9, 9, 11, 10, -1, -1, -1, -1, -1, -1, -1},
+        {9, 5, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {9, 5, 4, 0, 8, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 5, 4, 1, 5, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {8, 5, 4, 8, 3, 5, 3, 1, 5, -1, -1, -1, -1, -1, -1, -1},
+        {1, 2, 10, 9, 5, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {3, 0, 8, 1, 2, 10, 4, 9, 5, -1, -1, -1, -1, -1, -1, -1},
+        {5, 2, 10, 5, 4, 2, 4, 0, 2, -1, -1, -1, -1, -1, -1, -1},
+        {2, 10, 5, 3, 2, 5, 3, 5, 4, 3, 4, 8, -1, -1, -1, -1},
+        {9, 5, 4, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 11, 2, 0, 8, 11, 4, 9, 5, -1, -1, -1, -1, -1, -1, -1},
+        {0, 5, 4, 0, 1, 5, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1},
+        {2, 1, 5, 2, 5, 8, 2, 8, 11, 4, 8, 5, -1, -1, -1, -1},
+        {10, 3, 11, 10, 1, 3, 9, 5, 4, -1, -1, -1, -1, -1, -1, -1},
+        {4, 9, 5, 0, 8, 1, 8, 10, 1, 8, 11, 10, -1, -1, -1, -1},
+        {5, 4, 0, 5, 0, 11, 5, 11, 10, 11, 0, 3, -1, -1, -1, -1},
+        {5, 4, 8, 5, 8, 10, 10, 8, 11, -1, -1, -1, -1, -1, -1, -1},
+        {9, 7, 8, 5, 7, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {9, 3, 0, 9, 5, 3, 5, 7, 3, -1, -1, -1, -1, -1, -1, -1},
+        {0, 7, 8, 0, 1, 7, 1, 5, 7, -1, -1, -1, -1, -1, -1, -1},
+        {1, 5, 3, 3, 5, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {9, 7, 8, 9, 5, 7, 10, 1, 2, -1, -1, -1, -1, -1, -1, -1},
+        {10, 1, 2, 9, 5, 0, 5, 3, 0, 5, 7, 3, -1, -1, -1, -1},
+        {8, 0, 2, 8, 2, 5, 8, 5, 7, 10, 5, 2, -1, -1, -1, -1},
+        {2, 10, 5, 2, 5, 3, 3, 5, 7, -1, -1, -1, -1, -1, -1, -1},
+        {7, 9, 5, 7, 8, 9, 3, 11, 2, -1, -1, -1, -1, -1, -1, -1},
+        {9, 5, 7, 9, 7, 2, 9, 2, 0, 2, 7, 11, -1, -1, -1, -1},
+        {2, 3, 11, 0, 1, 8, 1, 7, 8, 1, 5, 7, -1, -1, -1, -1},
+        {11, 2, 1, 11, 1, 7, 7, 1, 5, -1, -1, -1, -1, -1, -1, -1},
+        {9, 5, 8, 8, 5, 7, 10, 1, 3, 10, 3, 11, -1, -1, -1, -1},
+        {5, 7, 0, 5, 0, 9, 7, 11, 0, 1, 0, 10, 11, 10, 0, -1},
+        {11, 10, 0, 11, 0, 3, 10, 5, 0, 8, 0, 7, 5, 7, 0, -1},
+        {11, 10, 5, 7, 11, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {10, 6, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 8, 3, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {9, 0, 1, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 8, 3, 1, 9, 8, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1},
+        {1, 6, 5, 2, 6, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 6, 5, 1, 2, 6, 3, 0, 8, -1, -1, -1, -1, -1, -1, -1},
+        {9, 6, 5, 9, 0, 6, 0, 2, 6, -1, -1, -1, -1, -1, -1, -1},
+        {5, 9, 8, 5, 8, 2, 5, 2, 6, 3, 2, 8, -1, -1, -1, -1},
+        {2, 3, 11, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {11, 0, 8, 11, 2, 0, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1},
+        {0, 1, 9, 2, 3, 11, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1},
+        {5, 10, 6, 1, 9, 2, 9, 11, 2, 9, 8, 11, -1, -1, -1, -1},
+        {6, 3, 11, 6, 5, 3, 5, 1, 3, -1, -1, -1, -1, -1, -1, -1},
+        {0, 8, 11, 0, 11, 5, 0, 5, 1, 5, 11, 6, -1, -1, -1, -1},
+        {3, 11, 6, 0, 3, 6, 0, 6, 5, 0, 5, 9, -1, -1, -1, -1},
+        {6, 5, 9, 6, 9, 11, 11, 9, 8, -1, -1, -1, -1, -1, -1, -1},
+        {5, 10, 6, 4, 7, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {4, 3, 0, 4, 7, 3, 6, 5, 10, -1, -1, -1, -1, -1, -1, -1},
+        {1, 9, 0, 5, 10, 6, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1},
+        {10, 6, 5, 1, 9, 7, 1, 7, 3, 7, 9, 4, -1, -1, -1, -1},
+        {6, 1, 2, 6, 5, 1, 4, 7, 8, -1, -1, -1, -1, -1, -1, -1},
+        {1, 2, 5, 5, 2, 6, 3, 0, 4, 3, 4, 7, -1, -1, -1, -1},
+        {8, 4, 7, 9, 0, 5, 0, 6, 5, 0, 2, 6, -1, -1, -1, -1},
+        {7, 3, 9, 7, 9, 4, 3, 2, 9, 5, 9, 6, 2, 6, 9, -1},
+        {3, 11, 2, 7, 8, 4, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1},
+        {5, 10, 6, 4, 7, 2, 4, 2, 0, 2, 7, 11, -1, -1, -1, -1},
+        {0, 1, 9, 4, 7, 8, 2, 3, 11, 5, 10, 6, -1, -1, -1, -1},
+        {9, 2, 1, 9, 11, 2, 9, 4, 11, 7, 11, 4, 5, 10, 6, -1},
+        {8, 4, 7, 3, 11, 5, 3, 5, 1, 5, 11, 6, -1, -1, -1, -1},
+        {5, 1, 11, 5, 11, 6, 1, 0, 11, 7, 11, 4, 0, 4, 11, -1},
+        {0, 5, 9, 0, 6, 5, 0, 3, 6, 11, 6, 3, 8, 4, 7, -1},
+        {6, 5, 9, 6, 9, 11, 4, 7, 9, 7, 11, 9, -1, -1, -1, -1},
+        {10, 4, 9, 6, 4, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {4, 10, 6, 4, 9, 10, 0, 8, 3, -1, -1, -1, -1, -1, -1, -1},
+        {10, 0, 1, 10, 6, 0, 6, 4, 0, -1, -1, -1, -1, -1, -1, -1},
+        {8, 3, 1, 8, 1, 6, 8, 6, 4, 6, 1, 10, -1, -1, -1, -1},
+        {1, 4, 9, 1, 2, 4, 2, 6, 4, -1, -1, -1, -1, -1, -1, -1},
+        {3, 0, 8, 1, 2, 9, 2, 4, 9, 2, 6, 4, -1, -1, -1, -1},
+        {0, 2, 4, 4, 2, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {8, 3, 2, 8, 2, 4, 4, 2, 6, -1, -1, -1, -1, -1, -1, -1},
+        {10, 4, 9, 10, 6, 4, 11, 2, 3, -1, -1, -1, -1, -1, -1, -1},
+        {0, 8, 2, 2, 8, 11, 4, 9, 10, 4, 10, 6, -1, -1, -1, -1},
+        {3, 11, 2, 0, 1, 6, 0, 6, 4, 6, 1, 10, -1, -1, -1, -1},
+        {6, 4, 1, 6, 1, 10, 4, 8, 1, 2, 1, 11, 8, 11, 1, -1},
+        {9, 6, 4, 9, 3, 6, 9, 1, 3, 11, 6, 3, -1, -1, -1, -1},
+        {8, 11, 1, 8, 1, 0, 11, 6, 1, 9, 1, 4, 6, 4, 1, -1},
+        {3, 11, 6, 3, 6, 0, 0, 6, 4, -1, -1, -1, -1, -1, -1, -1},
+        {6, 4, 8, 11, 6, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {7, 10, 6, 7, 8, 10, 8, 9, 10, -1, -1, -1, -1, -1, -1, -1},
+        {0, 7, 3, 0, 10, 7, 0, 9, 10, 6, 7, 10, -1, -1, -1, -1},
+        {10, 6, 7, 1, 10, 7, 1, 7, 8, 1, 8, 0, -1, -1, -1, -1},
+        {10, 6, 7, 10, 7, 1, 1, 7, 3, -1, -1, -1, -1, -1, -1, -1},
+        {1, 2, 6, 1, 6, 8, 1, 8, 9, 8, 6, 7, -1, -1, -1, -1},
+        {2, 6, 9, 2, 9, 1, 6, 7, 9, 0, 9, 3, 7, 3, 9, -1},
+        {7, 8, 0, 7, 0, 6, 6, 0, 2, -1, -1, -1, -1, -1, -1, -1},
+        {7, 3, 2, 6, 7, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {2, 3, 11, 10, 6, 8, 10, 8, 9, 8, 6, 7, -1, -1, -1, -1},
+        {2, 0, 7, 2, 7, 11, 0, 9, 7, 6, 7, 10, 9, 10, 7, -1},
+        {1, 8, 0, 1, 7, 8, 1, 10, 7, 6, 7, 10, 2, 3, 11, -1},
+        {11, 2, 1, 11, 1, 7, 10, 6, 1, 6, 7, 1, -1, -1, -1, -1},
+        {8, 9, 6, 8, 6, 7, 9, 1, 6, 11, 6, 3, 1, 3, 6, -1},
+        {0, 9, 1, 11, 6, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {7, 8, 0, 7, 0, 6, 3, 11, 0, 11, 6, 0, -1, -1, -1, -1},
+        {7, 11, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {7, 6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {3, 0, 8, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 1, 9, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {8, 1, 9, 8, 3, 1, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1},
+        {10, 1, 2, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 2, 10, 3, 0, 8, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1},
+        {2, 9, 0, 2, 10, 9, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1},
+        {6, 11, 7, 2, 10, 3, 10, 8, 3, 10, 9, 8, -1, -1, -1, -1},
+        {7, 2, 3, 6, 2, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {7, 0, 8, 7, 6, 0, 6, 2, 0, -1, -1, -1, -1, -1, -1, -1},
+        {2, 7, 6, 2, 3, 7, 0, 1, 9, -1, -1, -1, -1, -1, -1, -1},
+        {1, 6, 2, 1, 8, 6, 1, 9, 8, 8, 7, 6, -1, -1, -1, -1},
+        {10, 7, 6, 10, 1, 7, 1, 3, 7, -1, -1, -1, -1, -1, -1, -1},
+        {10, 7, 6, 1, 7, 10, 1, 8, 7, 1, 0, 8, -1, -1, -1, -1},
+        {0, 3, 7, 0, 7, 10, 0, 10, 9, 6, 10, 7, -1, -1, -1, -1},
+        {7, 6, 10, 7, 10, 8, 8, 10, 9, -1, -1, -1, -1, -1, -1, -1},
+        {6, 8, 4, 11, 8, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {3, 6, 11, 3, 0, 6, 0, 4, 6, -1, -1, -1, -1, -1, -1, -1},
+        {8, 6, 11, 8, 4, 6, 9, 0, 1, -1, -1, -1, -1, -1, -1, -1},
+        {9, 4, 6, 9, 6, 3, 9, 3, 1, 11, 3, 6, -1, -1, -1, -1},
+        {6, 8, 4, 6, 11, 8, 2, 10, 1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 2, 10, 3, 0, 11, 0, 6, 11, 0, 4, 6, -1, -1, -1, -1},
+        {4, 11, 8, 4, 6, 11, 0, 2, 9, 2, 10, 9, -1, -1, -1, -1},
+        {10, 9, 3, 10, 3, 2, 9, 4, 3, 11, 3, 6, 4, 6, 3, -1},
+        {8, 2, 3, 8, 4, 2, 4, 6, 2, -1, -1, -1, -1, -1, -1, -1},
+        {0, 4, 2, 4, 6, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 9, 0, 2, 3, 4, 2, 4, 6, 4, 3, 8, -1, -1, -1, -1},
+        {1, 9, 4, 1, 4, 2, 2, 4, 6, -1, -1, -1, -1, -1, -1, -1},
+        {8, 1, 3, 8, 6, 1, 8, 4, 6, 6, 10, 1, -1, -1, -1, -1},
+        {10, 1, 0, 10, 0, 6, 6, 0, 4, -1, -1, -1, -1, -1, -1, -1},
+        {4, 6, 3, 4, 3, 8, 6, 10, 3, 0, 3, 9, 10, 9, 3, -1},
+        {10, 9, 4, 6, 10, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {4, 9, 5, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 8, 3, 4, 9, 5, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1},
+        {5, 0, 1, 5, 4, 0, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1},
+        {11, 7, 6, 8, 3, 4, 3, 5, 4, 3, 1, 5, -1, -1, -1, -1},
+        {9, 5, 4, 10, 1, 2, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1},
+        {6, 11, 7, 1, 2, 10, 0, 8, 3, 4, 9, 5, -1, -1, -1, -1},
+        {7, 6, 11, 5, 4, 10, 4, 2, 10, 4, 0, 2, -1, -1, -1, -1},
+        {3, 4, 8, 3, 5, 4, 3, 2, 5, 10, 5, 2, 11, 7, 6, -1},
+        {7, 2, 3, 7, 6, 2, 5, 4, 9, -1, -1, -1, -1, -1, -1, -1},
+        {9, 5, 4, 0, 8, 6, 0, 6, 2, 6, 8, 7, -1, -1, -1, -1},
+        {3, 6, 2, 3, 7, 6, 1, 5, 0, 5, 4, 0, -1, -1, -1, -1},
+        {6, 2, 8, 6, 8, 7, 2, 1, 8, 4, 8, 5, 1, 5, 8, -1},
+        {9, 5, 4, 10, 1, 6, 1, 7, 6, 1, 3, 7, -1, -1, -1, -1},
+        {1, 6, 10, 1, 7, 6, 1, 0, 7, 8, 7, 0, 9, 5, 4, -1},
+        {4, 0, 10, 4, 10, 5, 0, 3, 10, 6, 10, 7, 3, 7, 10, -1},
+        {7, 6, 10, 7, 10, 8, 5, 4, 10, 4, 8, 10, -1, -1, -1, -1},
+        {6, 9, 5, 6, 11, 9, 11, 8, 9, -1, -1, -1, -1, -1, -1, -1},
+        {3, 6, 11, 0, 6, 3, 0, 5, 6, 0, 9, 5, -1, -1, -1, -1},
+        {0, 11, 8, 0, 5, 11, 0, 1, 5, 5, 6, 11, -1, -1, -1, -1},
+        {6, 11, 3, 6, 3, 5, 5, 3, 1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 2, 10, 9, 5, 11, 9, 11, 8, 11, 5, 6, -1, -1, -1, -1},
+        {0, 11, 3, 0, 6, 11, 0, 9, 6, 5, 6, 9, 1, 2, 10, -1},
+        {11, 8, 5, 11, 5, 6, 8, 0, 5, 10, 5, 2, 0, 2, 5, -1},
+        {6, 11, 3, 6, 3, 5, 2, 10, 3, 10, 5, 3, -1, -1, -1, -1},
+        {5, 8, 9, 5, 2, 8, 5, 6, 2, 3, 8, 2, -1, -1, -1, -1},
+        {9, 5, 6, 9, 6, 0, 0, 6, 2, -1, -1, -1, -1, -1, -1, -1},
+        {1, 5, 8, 1, 8, 0, 5, 6, 8, 3, 8, 2, 6, 2, 8, -1},
+        {1, 5, 6, 2, 1, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 3, 6, 1, 6, 10, 3, 8, 6, 5, 6, 9, 8, 9, 6, -1},
+        {10, 1, 0, 10, 0, 6, 9, 5, 0, 5, 6, 0, -1, -1, -1, -1},
+        {0, 3, 8, 5, 6, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {10, 5, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {11, 5, 10, 7, 5, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {11, 5, 10, 11, 7, 5, 8, 3, 0, -1, -1, -1, -1, -1, -1, -1},
+        {5, 11, 7, 5, 10, 11, 1, 9, 0, -1, -1, -1, -1, -1, -1, -1},
+        {10, 7, 5, 10, 11, 7, 9, 8, 1, 8, 3, 1, -1, -1, -1, -1},
+        {11, 1, 2, 11, 7, 1, 7, 5, 1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 8, 3, 1, 2, 7, 1, 7, 5, 7, 2, 11, -1, -1, -1, -1},
+        {9, 7, 5, 9, 2, 7, 9, 0, 2, 2, 11, 7, -1, -1, -1, -1},
+        {7, 5, 2, 7, 2, 11, 5, 9, 2, 3, 2, 8, 9, 8, 2, -1},
+        {2, 5, 10, 2, 3, 5, 3, 7, 5, -1, -1, -1, -1, -1, -1, -1},
+        {8, 2, 0, 8, 5, 2, 8, 7, 5, 10, 2, 5, -1, -1, -1, -1},
+        {9, 0, 1, 5, 10, 3, 5, 3, 7, 3, 10, 2, -1, -1, -1, -1},
+        {9, 8, 2, 9, 2, 1, 8, 7, 2, 10, 2, 5, 7, 5, 2, -1},
+        {1, 3, 5, 3, 7, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 8, 7, 0, 7, 1, 1, 7, 5, -1, -1, -1, -1, -1, -1, -1},
+        {9, 0, 3, 9, 3, 5, 5, 3, 7, -1, -1, -1, -1, -1, -1, -1},
+        {9, 8, 7, 5, 9, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {5, 8, 4, 5, 10, 8, 10, 11, 8, -1, -1, -1, -1, -1, -1, -1},
+        {5, 0, 4, 5, 11, 0, 5, 10, 11, 11, 3, 0, -1, -1, -1, -1},
+        {0, 1, 9, 8, 4, 10, 8, 10, 11, 10, 4, 5, -1, -1, -1, -1},
+        {10, 11, 4, 10, 4, 5, 11, 3, 4, 9, 4, 1, 3, 1, 4, -1},
+        {2, 5, 1, 2, 8, 5, 2, 11, 8, 4, 5, 8, -1, -1, -1, -1},
+        {0, 4, 11, 0, 11, 3, 4, 5, 11, 2, 11, 1, 5, 1, 11, -1},
+        {0, 2, 5, 0, 5, 9, 2, 11, 5, 4, 5, 8, 11, 8, 5, -1},
+        {9, 4, 5, 2, 11, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {2, 5, 10, 3, 5, 2, 3, 4, 5, 3, 8, 4, -1, -1, -1, -1},
+        {5, 10, 2, 5, 2, 4, 4, 2, 0, -1, -1, -1, -1, -1, -1, -1},
+        {3, 10, 2, 3, 5, 10, 3, 8, 5, 4, 5, 8, 0, 1, 9, -1},
+        {5, 10, 2, 5, 2, 4, 1, 9, 2, 9, 4, 2, -1, -1, -1, -1},
+        {8, 4, 5, 8, 5, 3, 3, 5, 1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 4, 5, 1, 0, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {8, 4, 5, 8, 5, 3, 9, 0, 5, 0, 3, 5, -1, -1, -1, -1},
+        {9, 4, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {4, 11, 7, 4, 9, 11, 9, 10, 11, -1, -1, -1, -1, -1, -1, -1},
+        {0, 8, 3, 4, 9, 7, 9, 11, 7, 9, 10, 11, -1, -1, -1, -1},
+        {1, 10, 11, 1, 11, 4, 1, 4, 0, 7, 4, 11, -1, -1, -1, -1},
+        {3, 1, 4, 3, 4, 8, 1, 10, 4, 7, 4, 11, 10, 11, 4, -1},
+        {4, 11, 7, 9, 11, 4, 9, 2, 11, 9, 1, 2, -1, -1, -1, -1},
+        {9, 7, 4, 9, 11, 7, 9, 1, 11, 2, 11, 1, 0, 8, 3, -1},
+        {11, 7, 4, 11, 4, 2, 2, 4, 0, -1, -1, -1, -1, -1, -1, -1},
+        {11, 7, 4, 11, 4, 2, 8, 3, 4, 3, 2, 4, -1, -1, -1, -1},
+        {2, 9, 10, 2, 7, 9, 2, 3, 7, 7, 4, 9, -1, -1, -1, -1},
+        {9, 10, 7, 9, 7, 4, 10, 2, 7, 8, 7, 0, 2, 0, 7, -1},
+        {3, 7, 10, 3, 10, 2, 7, 4, 10, 1, 10, 0, 4, 0, 10, -1},
+        {1, 10, 2, 8, 7, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {4, 9, 1, 4, 1, 7, 7, 1, 3, -1, -1, -1, -1, -1, -1, -1},
+        {4, 9, 1, 4, 1, 7, 0, 8, 1, 8, 7, 1, -1, -1, -1, -1},
+        {4, 0, 3, 7, 4, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {4, 8, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {9, 10, 8, 10, 11, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {3, 0, 9, 3, 9, 11, 11, 9, 10, -1, -1, -1, -1, -1, -1, -1},
+        {0, 1, 10, 0, 10, 8, 8, 10, 11, -1, -1, -1, -1, -1, -1, -1},
+        {3, 1, 10, 11, 3, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 2, 11, 1, 11, 9, 9, 11, 8, -1, -1, -1, -1, -1, -1, -1},
+        {3, 0, 9, 3, 9, 11, 1, 2, 9, 2, 11, 9, -1, -1, -1, -1},
+        {0, 2, 11, 8, 0, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {3, 2, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {2, 3, 8, 2, 8, 10, 10, 8, 9, -1, -1, -1, -1, -1, -1, -1},
+        {9, 10, 2, 0, 9, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {2, 3, 8, 2, 8, 10, 0, 1, 8, 1, 10, 8, -1, -1, -1, -1},
+        {1, 10, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {1, 3, 8, 9, 1, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 9, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {0, 3, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
+        {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}
+    };
+
+}
\ No newline at end of file