import java.applet.Applet;
import java.awt.*;
import java.awt.event.*;
//import java.util.Vector;
import java.util.Random;
import java.util.Arrays;
import java.util.Vector;
//import org.w3c.dom.Node.*;
//import org.w3c.dom.Document;
//import org.w3c.dom.Element;


public class eProb extends Panel
	implements Runnable, ActionListener {

  public static void main(String[] args) {
    int height = 400;
    int width = 400;
    Frame f = new Frame();
    f.addWindowListener(new java.awt.event.WindowAdapter() {
         public void windowClosing(java.awt.event.WindowEvent e) {
         System.exit(0);
         };
       });

    eProb panel = new eProb();
    f.add(panel);
    f.pack();
    panel.setSize(width, height); // same size as defined in the HTML APPLET
    panel.init();
    f.setSize(width, height + 20); // add 20, seems enough for the Frame title,
    f.setVisible(true);
  }

  final public static long serialVersionUID = 101;

  private int image_width, image_height;
  private double image_diameter;
  private Image zBuf;
  private Graphics offScreenGraphics;
  private Button go_1s_Button, go_2s_Button, go_2p_Button;
  private Button clearButton, pauseButton;
  private volatile Thread animationThread = null;
//  private Document doc = this.parent;
//  private Element nojavawarning = doc.getElementById("no_java_warning");
//  nojavawarning.setAttribute("foo", "bar");

  int blot_size = 1;
  private int bar_where = 20;
  private int bar_size = 2;     // thickness of scale bar / progress bar
  double kZoom = 20;		// pixels per Bohr radius
  Random generator = new Random();
  int blot_max = 10000;		// number of blots in scatter plot
  int blot_count = blot_max;
  String doingOrbital = "-none-";
  Color LtLtYellow = new Color(0xff, 0xff, 0xc0);
  Color LtLtBlue = new Color(0xef, 0xef, 0xff);
  int ntrpRails = 100;	// number of fence-rails in interpolation tables;
                        // number of fence-posts is one greater.
  double[] radius = new double[1+ntrpRails];
  double[] prob_r = new double[1+ntrpRails];
  double[] angle  = new double[1+ntrpRails];
  double[] prob_a = new double[1+ntrpRails];
  Vector<Button> frill = new Vector<Button>();

  public void init() {
    System.out.println("Starting....");
    image_width = getSize().width;
    image_height = getSize().height;
    image_diameter = Math.sqrt(Math.pow(image_width,2) + Math.pow(image_height,2));

    setBackground(Color.yellow);
    zBuf = createImage(image_width, image_height);
    offScreenGraphics = zBuf.getGraphics();
    // Automatic in some systems, not in others
    offScreenGraphics.setColor(Color.black);

    go_1s_Button = new Button("Go 1s");
    go_1s_Button.addActionListener(this);
    frill.add(go_1s_Button);
    add(go_1s_Button);

    go_2s_Button = new Button("Go 2s");
    go_2s_Button.addActionListener(this);
    frill.add(go_2s_Button);
    add(go_2s_Button);

    go_2p_Button = new Button("Go 2p");
    go_2p_Button.addActionListener(this);
    frill.add(go_2p_Button);
    add(go_2p_Button);

    clearButton = new Button("Clear");
    clearButton.addActionListener(this);
    add(clearButton);

    pauseButton = new Button("Pause");
    pauseButton.addActionListener(this);
    add(pauseButton);

    setCount(go_1s_Button, "1s");
    init_1s_table();
    do_go();
  }

// When the "start" button is pressed, start the
//  animation thread if it is not already started.
//
//  When the "pause" button is pressed, notify the
//  animation thread it should stop.

  void do_go(){
    pauseButton.setLabel("Pause");

    if (blot_count >= blot_max) {
// clear the screen and reset the counter
      clearImage();
      blot_count = 0;
    }
    if (animationThread == null || !animationThread.isAlive()) {
      animationThread = new Thread(this);
      animationThread.start();
    }
    repaint();
  }

  void clearImage(){
    offScreenGraphics.setColor(LtLtYellow);
    offScreenGraphics.fillRect(0,0,image_width,image_height);
// Scale bar / progress bar:
    offScreenGraphics.setColor(Color.red);
    offScreenGraphics.fillRect(bar_where, image_height-bar_where,
                (int)kZoom,bar_size);
  }


///////////////////////////////
  void init_1s_table(){
    init_table_common();
    for (int ii = 0; ii < radius.length; ii++) {
      double rr = radius[ii];
// The 1s orbital wavefunction is proportional to
//    exp(-r)
// The probability per unit volume is therefore
//    exp(-2*r)
// And we integrate that w.r.t rdr (not just dr) to get:
      prob_r[ii] = 1. - (2.*rr + 1.)*Math.exp(-2.*rr);
    }

    for (int ii = 0; ii < angle.length; ii++) {
      double aa = angle[ii];
      prob_a[ii] =  aa/2.0/Math.PI;
    }

    dump_tables();
  }

///////////////////////////////
  void init_2s_table(){
    init_table_common();
    for (int ii = 0; ii < radius.length; ii++) {
      double rr = radius[ii];
// The 2s orbital wavefunction is proportional to
//   (2-r) exp(-r/2)
// and the square of that is
//    (2-r)*(2-r)* exp(-r)
// so the integral (rdr) of the probability density is:
      prob_r[ii] = 1. + (-0.5*rr*rr*rr + 0.5*rr*rr - rr -1)*Math.exp(-rr);
    }

    for (int ii = 0; ii < angle.length; ii++) {
      double aa = angle[ii];
      prob_a[ii] =  aa/2.0/Math.PI;
    }

    dump_tables();
  }

////////////////////////////////////////
// the 2p orbital wavefunction is
//	 r exp(-r/2) cos(theta).
// the probability density is proportional to
//       r^2 exp(-r) cos^2 (theta)
// and the integral (rdr) of the radial piece is
// (rr*rr*rr/-6. - 0.5*rr*rr - rr - 1)*exp(-rr)
  void init_2p_table(){
    init_table_common();
    for (int ii = 0; ii < radius.length; ii++) {
      double rr = radius[ii];
      prob_r[ii] = 1. + (rr*rr*rr/-6. - 0.5*rr*rr - rr - 1)*Math.exp(-rr);
    }

    for (int ii = 0; ii < angle.length; ii++) {
      double aa = angle[ii];
      prob_a[ii] = (.5*Math.sin(2.*aa) + aa)/2.0/Math.PI;
    }

    dump_tables();
  }

  void init_table_common(){
    double aStep = 2.*Math.PI / ntrpRails;
    double rStep = image_diameter / kZoom / ntrpRails;
    for (int ii = 0; ii <= ntrpRails; ii++) {
      angle[ii]  = aStep * ii;
      radius[ii] = rStep * ii;
    }
  }

  void dump_tables(){
    boolean verbose = false;
    if (verbose) {
      System.out.print("aa0: "); System.out.println(angle[0]);
      System.out.print("pa0: "); System.out.println(prob_a[0]);
      System.out.print("aa-: "); System.out.println(angle[angle.length-2]);
      System.out.print("pa-: "); System.out.println(prob_a[angle.length-2]);
      System.out.print("aa+: "); System.out.println(angle[angle.length-1]);
      System.out.print("pa+: "); System.out.println(prob_a[angle.length-1]);


      System.out.print("rr0: "); System.out.println(radius[0]);
      System.out.print("pr0: "); System.out.println(prob_r[0]);
      System.out.print("rr-: "); System.out.println(radius[radius.length-2]);
      System.out.print("pr-: "); System.out.println(prob_r[radius.length-2]);
      System.out.print("rr+: "); System.out.println(radius[radius.length-1]);
      System.out.print("pr+: "); System.out.println(prob_r[radius.length-1]);
    }
  }

  void setCount(Button B, String shouldBe){
    for (int ii=0; ii < frill.size(); ii++) {
      // show all other buttons as inactive:
      Button other = frill.elementAt(ii);
      other.setBackground(Color.yellow);
    }
    // show "this" button as active:
    B.setBackground(LtLtBlue);
    if (!doingOrbital.equals(shouldBe)) {
      blot_count = blot_max;    // force a clearscreen
      doingOrbital = shouldBe;
    }
  }
  public void actionPerformed(ActionEvent event) {
    Button B = (Button)event.getSource();
    if (B == go_1s_Button) {
      setCount(B, "1s");
      init_1s_table();
      do_go();
    } else if (B == go_2s_Button) {
      setCount(B, "2s");
      init_2s_table();
      do_go();
    } else if (B == go_2p_Button) {
      setCount(B, "2p");
      init_2p_table();
      do_go();
    } else if (B == clearButton) {
// If animation thread already running, cause it
// to start over and keep running.
// If stopped, clear the screen and leave it stopped.
      clearImage();
      blot_count = 0;
      repaint();
    } else if (B == pauseButton) {
      animationThread = null;
      repaint();
    } else {
      System.out.print("eProb: Oops, unimplemented button: ");
      System.out.println(B);
    }
  }

// This only runs as part of the animation thread,
// not the parent thread.
// Runs until thread is killed.

  public void run() {

//    offScreenGraphics.drawString("Started.", 10, 10);
    Thread goodThread = Thread.currentThread();

    while(true) {
      if (blot_count >= blot_max) {
        // Changing the button label is cute,
	// but too subtle.
	//xx  pauseButton.setLabel("Done");
	break;
      }
      if (animationThread != goodThread) break;
      offScreenGraphics.setColor(Color.black);
      double xx, yy, rr, aa;
      double pr = generator.nextDouble();
      double pa = generator.nextDouble();
      try {
	rr = interpolate(pr, prob_r, radius);
      } catch (AssertionError foo) {
	System.out.print("%radius would extrapolate : pr = ");
	System.out.println(pr);
        rr = -1;
//xx	stop();
      }
      try {
        aa = interpolate(pa, prob_a, angle);
      } catch (AssertionError foo) {
	System.out.print("?angle would extrapolate! : pa = ");
	System.out.println(pa);
	rr = -1;
	aa = 0;
      }
      if (rr >= 0.){
	xx = rr * Math.cos(aa);
	yy = rr * Math.sin(aa);

	blot(xx,yy);
	offScreenGraphics.fillRect(bar_where, image_height-bar_where,
	  (int)(.5 + (kZoom*blot_count)/blot_max), bar_size);
	repaint();
	blot_count++;
      }
      pause(5);
    }
//    offScreenGraphics.drawString("Done.", 50, 10);
    repaint();
//xx    stop();
  }

// Call with xx and yy in atomic units,
// i.e. units where the Bohr radius is 1.0.
  void blot(double xx, double yy) {
    xx = xx*kZoom + image_width/2;
    yy = yy*kZoom + image_height/2;
    offScreenGraphics.fillRect((int)xx, (int)yy, blot_size,blot_size);
  }

// Override default update function.
// Skip the usual screen-clearing step of update
//  so that there is no "flicker" between each
//  drawing step.
  public void update(Graphics g) {
    paint(g);
  }


// Draw the  offscreen pixmap onto the applet window.
  public void paint(Graphics g) {
    g.drawImage(zBuf, 0, 0, this);
  }


// Sleep for the specified amount of time.
  private void pause(int milliseconds) {
    try {
      Thread.sleep(milliseconds);
    }
    catch(InterruptedException ie) {
      /* do nothing */
    }
  }

//////////////////
// Interpolation stuff.
//
// throws AssertionError if an extrapolation would be required
//
  double interpolate(double xx, double[] abs, double[] ord)
    throws AssertionError
  {
    int ii = Arrays.binarySearch(abs, xx);
    boolean verbose = false;
    if (ii >= 0) return ord[ii];
    ii = -(1+ii);	// now insertion point
    if (verbose && ii>=abs.length || ii>=ord.length) {
      System.out.print("xx: "); System.out.println(xx);
      System.out.print("a0: "); System.out.println(abs[0]);
      System.out.print("a+: "); System.out.println(abs[abs.length-1]);
      System.out.print("b0: "); System.out.println(ord[0]);
      System.out.print("b+: "); System.out.println(ord[abs.length-1]);
      System.out.print("ii: "); System.out.println(ii);
    }
    if (ii>=abs.length || ii>=ord.length) {
      throw new AssertionError();
    }
    return interp(xx, abs[ii-1], abs[ii], ord[ii-1], ord[ii]);
  }

  double interp(double xx, double aa, double AA, double bb, double BB){
    return ((xx-aa) * BB + (AA-xx) * bb) / (AA - aa);
  }

}

// Interpolation in java:
// http://java.sun.com/javase/6/docs/api/java/util/Arrays.html

// Solutions to the Schroedinger equation
// http://galileo.phys.virginia.edu/classes/751.mf1i.fall02/HydrogenAtom.htm

// Java Console:
// https://addons.mozilla.org/en-US/firefox/addon/141

// Maxima manual
// http://www.ma.utexas.edu/maxima/

// Maxima integration
// http://www.ma.utexas.edu/maxima/maxima_19.html

// Converting applet to stand-alone
// http://www.rgagnon.com/javadetails/java-0305.html