import java.awt.*;
import java.applet.*;
import java.awt.image.*;

/**
 * Test applet for using a triangular underlying mesh instead
 * of a square one
 */
public class Trid extends Applet
{
   final int COARSE_MESH_DIM = 1;
   final int MAX_LEVELS = 9;
   final int SIDE_DIM = COARSE_MESH_DIM * ( 1 << MAX_LEVELS );
   final double alpha = 1.0;

   public boolean generateDone = false;

   public boolean useBadInterpolation = false;
   // public boolean useBadInterpolation = true;

   public boolean debugTiming = true;

   public boolean useMIS = true;

   public final static int SHADE_HEIGHTMAP = 0;

   public final static int SHADE_SURFACE_NORMAL = 1;

   public int shadingMode = SHADE_SURFACE_NORMAL;
   // public int shadingMode = SHADE_HEIGHTMAP;

   public float [][] z;
   public float [][] zshade;
   public float [][][] norm;
   public boolean [][] zset;
   public int [][] iz;

   public void init()
   {
      generateHierarchicalFunction();
      Thread animator = new Thread() {
         public void run()
         {
            while ( true ) {
               animateCycle();
               try {
                  Thread.currentThread().sleep( 100 );
               } catch ( InterruptedException ie ) {
               }
            }
         }
      };
      animator.start();
   }

   public float randomFloat()
   {
      return ( float ) Math.random();
   }

   static float [] ztmp = new float[ 7 ];

   public void computeSurfaceNormal( int i, int j )
   {
      // fill in the <Ztmp> array with the surrounding points if
      // available or z[i][j] if not

      boolean allZSet = true;
      for ( int k = 1; k < 7; k++ ) ztmp[ k ] = z[ i ][ j ];
      if ( j > 0 ) {
         ztmp[ 1 ] = z[ i ][ j - 1 ];
         if ( !zset[ i ][ j - 1 ] ) allZSet = false;
      }
      if ( j < SIDE_DIM - i ) {
         ztmp[ 4 ] = z[ i ][ j + 1 ];
         if ( !zset[ i ][ j + 1 ] ) allZSet = false;
      }
      if ( i > 0 ) {
         if ( j > 0 ) {
            ztmp[ 6 ] = z[ i - 1 ][ j - 1 ];
            if ( !zset[ i - 1 ][ j - 1 ] ) allZSet = false;
         }
         if ( j < SIDE_DIM - ( i - 1 ) ) {
            ztmp[ 5 ] = z[ i - 1 ][ j ];
            if ( !zset[ i - 1 ][ j ] ) allZSet = false;
         }
      }
      if ( i < SIDE_DIM ) {
         if ( j + 1 < SIDE_DIM - ( i + 1 ) ) {
            ztmp[ 3 ] = z[ i + 1 ][ j + 1 ];
            if ( !zset[ i + 1 ][ j + 1 ] ) allZSet = false;
         }
         if ( j < SIDE_DIM - ( i + 1 ) ) {
            ztmp[ 2 ] = z[ i + 1 ][ j ];
            if ( !zset[ i + 1 ][ j ] ) allZSet = false;
         }
      }

      float zavg = ( ztmp[ 1 ] + ztmp[ 2 ] + ztmp[ 3 ] + 
         ztmp[ 4 ] + ztmp[ 5 ] + ztmp[ 6 ] ) / 6;
      for ( int k = 1; k <= 6; k++ ) ztmp[ k ] -= zavg;

      float bb = ( 2.0f * ( ztmp[ 4 ] - ztmp[ 1 ] ) + ztmp[ 3 ] + ztmp[ 5 ] -
         ztmp[ 2 ] - ztmp[ 6 ] ) / 8;
      float cc = ( ztmp[ 2 ] + ztmp[ 3 ] - ztmp[ 5 ] - ztmp[ 6 ] ) / 4;

      float b = 2.0f * bb;
      float c = sq3o2i * cc;
      float normalize = ( float ) Math.sqrt( 1.0 + b * b + c * c );

      norm[ i ][ j ][ 0 ] = - c / normalize;
      norm[ i ][ j ][ 1 ] = - b / normalize;
      norm[ i ][ j ][ 2 ] = 1.0f / normalize;

      if ( !allZSet ) {
         norm[ i ][ j ][ 0 ] = 0;
         norm[ i ][ j ][ 1 ] = 0;
         norm[ i ][ j ][ 2 ] = 0;
      }
   }

   public float computeSurfaceNormalDotLight( int i, int j ) {
      return norm[ i ][ j ][ 0 ] * xLight + norm[ i ][ j ][ 1 ] * yLight 
         + norm[ i ][ j ][ 2 ] * zLight;
   }

   final static float sq3o2i = ( float ) ( 2.0 / Math.sqrt( 3.0 ) );

   float xLightAngle = ( float ) ( 0.4 * Math.PI );
   float yLightAngle = ( float ) ( 0.3 * Math.PI );
   float xLight = ( float ) ( Math.cos( xLightAngle ) * Math.sin( yLightAngle ) );
   float yLight = ( float ) ( Math.sin( xLightAngle ) * Math.sin( yLightAngle ) );
   float zLight = ( float ) Math.cos( yLightAngle );

   // interp_weight is the relative weight of the nearer vertices 
   // compared to the farther vertices.  There are several mathematically
   // justifiable choices for interp_weight:
   //
   // Inf     this is equivalent to useBadInterpolation = true; the far
   //         vertices aren't weighted at all, which produces really
   //         extreme mesh artifacts
   //
   // 0       this uses ONLY the far vertices, and while it produces almost
   //         no grid artifacts it makes a REALLY weird picture.  The 
   //         height map is really gentle.
   //
   // 1       if we do the interpolation by least-squares fit to a plane,
   //         z(x,y) = A + B x + C y, it produces this choice.  It doesn't
   //         produce a mesh-orientation-independent picture, especially
   //         for alpha >= 2, and I don't know why.  Seems like it should.
   //
   // 9       if instead of a plane we fit a biquadratic with no x-y cross
   //         term (which means grid-orientation-dependence!) we get this
   //         weight.  It works ok, but there's definite grid dependence at
   //         large values of alpha
   //
   // sqrt(3) this one isn't mathematically justified but it seems to work
   //         pretty well.  2 and 3 also work reasonably well, producing only
   //         the "nipples" effect at alpha >~= 2.
   //
   // Is the nipple effect necessarily bad?  If the soil has a constant
   // angle of repose, then alpha should really be around 1.  It's really
   // more complicated than this; I need to think about it.  alpha = 2 is
   // probably way too big.  A realistic calculation would be alpha = 1 
   // with some sort of "maximum angle of repose" smoothing happening 
   // at each step.

   float interp_weight = 2.0f;
   float interp_weightx2 = interp_weight * 2;

   public void generateHierarchicalFunction()
   {
      // float base_height = 1.0f;
      // float height = base_height / 2;
      long start_time = System.currentTimeMillis();
      float base_height = 0;
      float height = 1.0f;

      double base_angle = Math.PI / 6.0;
      height = ( float ) ( Math.tan( base_angle ) * ( 1 << MAX_LEVELS ) );

      // fill in the top-level random function on
      // the triangular coarse mesh

      int stride = 1 << MAX_LEVELS;
      z = new float[ SIDE_DIM + 1 ][];
      norm = new float[ SIDE_DIM + 1 ][][];
      iz = new int[ SIDE_DIM + 1][];
      zset = new boolean[ SIDE_DIM + 1 ][];
      if ( shadingMode == SHADE_SURFACE_NORMAL ) {
         zshade = new float[ SIDE_DIM + 1 ][];
      }
      for ( int i = 0; i <= SIDE_DIM; i++ ) {
         z[ i ] = new float[ SIDE_DIM - i + 1 ];
         norm[ i ] = new float[ SIDE_DIM - i + 1 ][ 3 ];
         zset[ i ] = new boolean[ SIDE_DIM - i + 1 ];
         iz[ i ] = new int[ SIDE_DIM - i + 1 ];
         if ( shadingMode == SHADE_SURFACE_NORMAL ) {
            zshade[ i ] = new float[ SIDE_DIM - i + 1 ];
         }
         for ( int j = 0; j <= SIDE_DIM - i; j++ ) {
            zset[ i ][ j ] = false;
         }
      }
      for ( int i = 0; i <= SIDE_DIM; i += stride ) {
         for ( int j = 0; j <= SIDE_DIM - i; j += stride ) {
            z[ i ][ j ] = base_height * ( randomFloat() - 0.5f );
            zset[ i ][ j ] = true;
         }
      }

      // perform MAX_LEVELS passes of interpolate / perturb over the mesh

      for ( int l = MAX_LEVELS; l > 0; l-- ) {
         stride = 1 << l;
         int half_stride = stride / 2;
         for ( int i = half_stride; i <= SIDE_DIM; i += stride ) {
            for ( int j = 0; j <= SIDE_DIM; j += stride ) {

               // pass 1 interp: x = j, y = i, interp from (x, y +- half)
               // (Good: interp from (x +- 1, y) and (x, y +- half) )

               if ( i + half_stride <= SIDE_DIM - j ) {
                  z[ j ][ i ] = ( z[ j ][ i - half_stride ] +
                     z[ j ][ i + half_stride ] ) * interp_weight;
                  float npts = interp_weightx2;
                  if ( !useBadInterpolation ) {
                     if ( ( j > 0 ) && ( i + half_stride <=
                        SIDE_DIM - ( j - stride ) ) ) {
                        z[ j ][ i ] += z[ j - stride ][ i + half_stride ];
                        npts++;
                     }
                     if ( ( j + stride <= SIDE_DIM ) && ( i - half_stride <=
                        SIDE_DIM - ( j + stride ) ) ) {
                        z[ j ][ i ] += z[ j + stride ][ i - half_stride ];
                        npts++;
                     }
                  }
                  z[ j ][ i ] /= npts;
                  z[ j ][ i ] += height * ( randomFloat() - 0.5f );
                  zset[ j ][ i ] = true;
               }

               // pass 2 interp: y = j, x = i, interp from (x +- half, y)

               if ( j <= SIDE_DIM - ( i + half_stride ) ) {
                  z[ i ][ j ] = ( z[ i - half_stride ][ j ] +
                     z[ i + half_stride ][ j ] ) * interp_weight;
                  float npts = interp_weightx2;
                  if ( !useBadInterpolation ) {
                     if ( ( i + half_stride <= SIDE_DIM ) && ( i - half_stride 
                        + j <= SIDE_DIM ) && ( j >= stride ) ) {
                        z[ i ][ j ] += z[ i + half_stride ][ j - stride ];
                        npts++;
                     }
                     if ( ( i > 0 ) && ( i + j + half_stride <= SIDE_DIM ) ) {
                        z[ i ][ j ] += z[ i - half_stride ][ j + stride ];
                        npts++;
                     }
                  }
                  z[ i ][ j ] /= npts;
                  z[ i ][ j ] += height * ( randomFloat() - 0.5f );
                  zset[ i ][ j ] = true;
               }

               // pass 3 interp: z = j, x = i, interp from (x +- half, z)
               //   where y = DIM - x - z

               int k = SIDE_DIM - i - j;
               int k1 = k - half_stride;
               int k2 = k + half_stride;
               if ( ( k1 >= 0 ) && ( k2 <= SIDE_DIM - i ) ) {
                  z[ i ][ k ] = ( z[ i - half_stride ][ k2 ] +
                     z[ i + half_stride ][ k1 ] ) * interp_weight;
                  float npts = interp_weightx2;
                  if ( !useBadInterpolation ) {
                     if ( ( i >= half_stride ) && ( i - half_stride + k1 <=
                        SIDE_DIM ) ) {
                        z[ i ][ k ] += z[ i - half_stride ][ k1 ];
                        npts++;
                     }
                     if ( ( i + half_stride <= SIDE_DIM ) && ( i + half_stride
                        + k2 <= SIDE_DIM ) && ( k2 >= 0 ) ) {
                        z[ i ][ k ] += z[ i + half_stride ][ k2 ];
                        npts++;
                     }
                  }
                  z[ i ][ k ] /= npts;
                  z[ i ][ k ] += height * ( randomFloat() - 0.5f );
                  zset[ i ][ k ] = true;
               }
            }
         }
         height *= ( float ) Math.pow( 2.0, -alpha );
      }

      long hf_end_time = System.currentTimeMillis();

      if ( shadingMode == SHADE_SURFACE_NORMAL ) {
         for ( int i = 0; i <= SIDE_DIM; i++ ) {
            for ( int j = 0; j <= SIDE_DIM - i; j++ ) {
               computeSurfaceNormal( i, j );
            }
         }
      }
      long normal_end_time = System.currentTimeMillis();
      if ( debugTiming ) {
         System.out.println( "HF Gen took " + ( 0.001 * ( hf_end_time -
            start_time ) ) + " s" );
         System.out.println( "Normals took " + ( 0.001 * ( normal_end_time
            - hf_end_time ) ) + " s" );
      }
      generateDone = true;
      repaint();
   }

   public void reshade()
   {
      if ( shadingMode == SHADE_SURFACE_NORMAL ) {
         for ( int i = 0; i <= SIDE_DIM; i++ ) {
            for ( int j = 0; j <= SIDE_DIM - i; j++ ) {
               zshade[ i ][ j ] = computeSurfaceNormalDotLight( i, j );
               if ( zshade[ i ][ j ] < 0 ) {
                  zshade[ i ][ j ] = 0;
               }
            }
         }
      }
   }

   double theta = 0.0;

   public void animateCycle()
   {
      theta += 0.02 * Math.PI;
      if ( theta > 2.0 * Math.PI ) theta -= 2.0 * Math.PI;
      xLightAngle = ( float ) Math.cos( theta );
      yLightAngle = ( float ) Math.sin( theta );
/*
      xLightAngle += ( float ) ( 0.02 * Math.PI );
      if ( xLightAngle > Math.PI - 0.1 ) {
         xLightAngle = 0.1f;
      }
      yLightAngle += ( float ) ( 0.01 * Math.PI );
      if ( yLightAngle > Math.PI - 0.1 ) {
         yLightAngle = 0.1f;
      }
*/
      xLight = ( float ) xLightAngle;
      yLight = ( float ) yLightAngle;
      zLight = ( float ) 1.0f;
      float norm = ( float ) Math.sqrt( xLight * xLight + yLight * yLight +
         1.0 );
      xLight /= norm;
      yLight /= norm;
      zLight /= norm;
      reshade();
      repaint();
   }

   boolean colorsDone = false;

   Color [] colors = new Color[ 100 ];

   float zmin, zmax;

   Image buffer = null;

   int [] pix = new int[ ( SIDE_DIM + 1 ) * ( SIDE_DIM + 1 ) ];

   boolean repaintInProgress = false;

   MemoryImageSource mis = null;

   Image misbuf = null;

   public void paint( Graphics gg )
   {
      if ( repaintInProgress ) {
         return;
      }
      repaintInProgress = true;
      if ( generateDone ) {
         if ( !colorsDone ) {
            long color_start_time = System.currentTimeMillis();
            zmax = -1000;
            zmin = 1000;
            for ( int i = 0; i <= SIDE_DIM; i++ ) {
               for ( int j = 0; j <= SIDE_DIM - i; j++ ) {
                  if ( shadingMode == SHADE_SURFACE_NORMAL ) {
                     if ( zshade[ i ][ j ] > zmax ) {
                        zmax = zshade[ i ][ j ];
                     } 
                     if ( zshade[ i ][ j ] < zmin ) {
                        zmin = zshade[ i ][ j ];
                     }
                  } else {
                     if ( z[ i ][ j ] > zmax ) {
                        zmax = z[ i ][ j ];
                     } 
                     if ( z[ i ][ j ] < zmin ) {
                        zmin = z[ i ][ j ];
                     }
                  }
               }
            }

            for ( int i = 0; i < 100; i++ ) {
               colors[ i ] = new Color( 2 * i, 2 * i, 2 * i );
            }
            long color_end_time = System.currentTimeMillis();
            if ( debugTiming ) {
               System.out.println( "Color gen took " + ( 0.001 * (
                  color_end_time - color_start_time ) ) + " s" );
            }
            colorsDone = true;
         }

         if ( buffer == null ) {
            Dimension dim = getSize();
            buffer = createImage( dim.width, dim.height );
         }
         Graphics g = buffer.getGraphics();

         // draw it

         int scale = 1 + ( 400 / ( SIDE_DIM + 1 ) );
         // int ypad = 100;
         int ypad = 10;
         int eff_height = ( SIDE_DIM + 1 ) * scale + ypad;

         // g.setColor( new Color( 128, 128, 255 ) );
         g.setColor( Color.black );
         g.fillRect( 0, 0, eff_height, eff_height );

         long draw_start_time = System.currentTimeMillis();
         float ndz = 99 / ( zmax - zmin );
         
         int bluepix = ( 255 << 24 ) + ( 128 << 16 ) + ( 128 << 8 ) +
            255;
         bluepix = 0;
         for ( int i = 0; i < pix.length; i++ ) pix[ i ] = bluepix;
         
         int ct = 0;
         int mask = ( 1 << 16 ) + ( 1 << 8 ) + 1;
         for ( int i = 0; i <= SIDE_DIM; i++ ) {
            ct += i / 2;
            for ( int j = 0; j <= SIDE_DIM - i; j++ ) {
               if ( shadingMode == SHADE_SURFACE_NORMAL ) {
                  iz[ i ][ j ] = ( int ) ( ndz * ( zshade[ i ][ j ] - zmin ) );
               } else {
                  iz[ i ][ j ] = ( int ) ( ndz * ( z[ i ][ j ] - zmin ) );
               }
               pix[ ct++ ] = ( 2 * iz[ i ][ j ] * mask ) + ( 255 << 24 );
            }
            ct += i - i / 2;
         }
         long draw_mid_time = System.currentTimeMillis();

         if ( ( useMIS ) && ( scale == 1 ) ) {
            if ( mis == null ) {
               mis = new MemoryImageSource( SIDE_DIM + 1,
                  SIDE_DIM + 1, pix, 0, SIDE_DIM + 1 );
               mis.setAnimated( true );
               misbuf = createImage( mis );
            } else {
               mis.newPixels( 0, 0, SIDE_DIM + 1, SIDE_DIM + 1 );
            }
            g.drawImage( misbuf, 0, ypad, null );
         } else {
            for ( int c = 0; c < 100; c++ ) {
               g.setColor( colors[ c ] );
               int ys = eff_height;
               for ( int i = 0; i <= SIDE_DIM; i++ ) {
                  for ( int j = 0; j <= SIDE_DIM - i; j++ ) {
                     if ( ( zset[ i ][ j ] ) && ( iz[ i ][ j ] == c ) ) {
                        int xs = scale * ( j + i / 2 );
                        g.fillRect( xs, ys - scale, scale, scale );
                     }
                  }
                  ys -= scale;
               }
            }
         }
         long draw_end_time = System.currentTimeMillis();
         if ( debugTiming ) {
            System.out.println( "iz calc took " + ( 0.001 * ( draw_mid_time -
               draw_start_time ) ) + " s" );
            System.out.println( "Draw took " + ( 0.001 * ( draw_end_time -
               draw_mid_time ) ) + " s" );
            System.out.println( "   scale = " + scale );
         }

         gg.drawImage( buffer, 0, 0, null );
      }
      repaintInProgress = false;
   }

   public void update( Graphics g ) {
      paint( g );
   }
}