import ccj.*;

public class Fractal extends GraphicsApplet
{  public void run()
   {  double[][][] transforms = new double[4][3][3];

      /* set up transforms */
      transforms[0][0][0] = 0.05;
      transforms[0][0][1] = 0;
      transforms[0][0][2] = 0;
      transforms[0][1][0] = 0;
      transforms[0][1][1] = 0.16;
      transforms[0][1][2] = -4.2;

      transforms[1][0][0] = 0.85;
      transforms[1][0][1] = 0.04;
      transforms[1][0][2] = 0.20;
      transforms[1][1][0] = -0.04;
      transforms[1][1][1] = 0.85;
      transforms[1][1][2] = 0.85;

      transforms[2][0][0] = 0.2;
      transforms[2][0][1] = -0.26;
      transforms[2][0][2] = -1.3;
      transforms[2][1][0] = 0.23;
      transforms[2][1][1] = 0.22;
      transforms[2][1][2] = -2.3;

      transforms[3][0][0] = -0.15;
      transforms[3][0][1] = 0.28;
      transforms[3][0][2] = 1.4;
      transforms[3][1][0] = 0.26;
      transforms[3][1][1] = 0.24;
      transforms[3][1][2] = -3.36;

      /* fill the third row of every transform with 0 0 1 */
      int i;
      for (i = 0; i < transforms.length; i++)
      {  transforms[i][2][0] = 0;
         transforms[i][2][1] = 0;
         transforms[i][2][2] = 1;
      }

      /* compute probabilities */
      double[] prob = new double[transforms.length];
                /* probability to select the ith transformation */
      double sum = 0;
      for (i = 0; i < transforms.length; i++)
      {  prob[i] = Math.abs(determinant(transforms[i]));
         sum = sum + prob[i];
      }
      for (i = 0; i < transforms.length; i++)
         prob[i] = prob[i] / sum;

      double[] p = new double[3];
         /* the point to be transformed */
      p[0] = 1; p[1] = 1; p[2] = 1;

      /* compute and plot transforms of p */
      int t;
      for (t = 1; t <= NTRIES; t++)
      {  double r = Numeric.randomDouble(0, 1);
         i = 0;
         while (r > prob[i] && i < transforms.length - 1)
         {  r = r - prob[i];
            i++;
         }
         /* now pick the ith transform */
         p = vproduct(transforms[i], p);
         new Point(p[0], p[1]).draw();
      }
      
   }

   public static double[] vproduct(double[][] a, double[] v)
   {  double[] r = new double[a.length];

      int i;
      for (i = 0; i < a.length; i++)
      {  if (a[i].length != v.length)
            throw new IllegalArgumentException();
         int k;
         double sum = 0;
         for (k = 0; k < v.length; k++)
            sum = sum + a[i][k] * v[k];
         r[i] = sum;
      }
      return r;
   }

   public static double determinant(double[][] m)
   {  if (m.length != 3 || m[0].length != 3
            || m[2][0] != 0 || m[2][1] != 0 || m[2][2] != 1)
         throw new IllegalArgumentException();         
      return m[0][0] * m[1][1] - m[0][1] * m[1][0];
   }

   private static final int NTRIES = 800;
}