/**
 * width 401, height 301
 * Rutherford scattering from a fixed force center for impact parameters
 * less than five times distance of closest approach at head-on impact
 * SEQ - 26 particles at regular impact parameter increments
 * RND - 26 particles randomly distributed across the area of the beam
 * CNT - angular distribution for 500 randomly distributed particles
 * 
 * The incident trajectory is specified by two parameters: a, the distance of 
 * closest approach in a head-on collision, and b, the impact parameter. In 
 * the example, the incident beam has a circular cross-section 
 * of radius b = 2.5a. The starting point for the numerical calculation is 40a 
 * from the scattering center (twice the distance shown in the plot). The program 
 * offers three choices: (a) trajectories for impact parameters stepped uniformly 
 * across the beam, (b) trajectories for the same number of incident particles 
 * distributed randomly across the area of the beam, and (c) a histogram versus 
 * scattering angle for 500 randomly distributed incident particles. (The 
 * histogram scale is multiplied by 10 for angles greater than 60 degrees. This 
 * multiplication distorts the visual interpretation, but the total number of 
 * counts in each scale region is shown on screen.) The histogram also shows the 
 * distribution calculated by analytic methods. The analytic treatment evaluates 
 * the integral in the Apsidal Angles section and finds that the impact parameter 
 * and scattering angle are related by tan(theta/2) = a/2b. (We will use z for 
 * this expression - an expression that predicts, for example, that no particles 
 * in our beam will be scattered through less than 22.6 degrees.) The angular 
 * distribution calculated using this expression, the Rutherford scattering 
 * formula, is usually written in terms of the inverse fourth power of a sine 
 * function, but when all the particles are being counted it is more convenient 
 * to write the angular distribution as proportional to (1+z^2)/z^3. This is the 
 * function plotted on the histogram. It predicts that 60 of the 500 particles 
 * will be scattered through more than 60 degrees. 
 * 
 * HTML file:   <applet code= fscat.class width=401 height= 301></applet>
 * 
 * C.A. Bertulani, 09/05/2000
**/

import java.applet.Applet;
import java.awt.*;
import java.awt.event.ActionListener;
import java.awt.event.ActionEvent;
public class fscat extends Applet 
                  implements ActionListener {
  // stepped impact parameter, kk =1, random impact parameter, kk=2,
  // and counter of scattered particles, kk=3.	
  int kk=0;  
  // Declarations of Buttons
  String b1s="SEQ"; Button b1=new Button(b1s);
  String b2s="RND"; Button b2=new Button(b2s);
  String b3s="CNT"; Button b3=new Button(b3s);
  public void init() {
	// FlowLayout puts components in a row, sized at their preferred size. 
	// If the horizontal space in the container is too small to put all the 
	// components in one row, FlowLayout uses multiple rows. Within each row, 
	// components are left-aligned when the FlowLayout is created. 
    setLayout(new FlowLayout(FlowLayout.LEFT));
	// add buttons
    add(b1); add(b2); add(b3);
	// listen to action on buttons
    b1.addActionListener(this);
    b2.addActionListener(this);
    b3.addActionListener(this);
    setBackground(Color.white);
  }
  double hit(double y, Graphics g) {        
	// Feynman Algorithm with trajectory
    int phase, x1, y1, x2, y2;                         
    double x, vx, vy, dt, r3;
	// initial conditions for positions and velocities 
    x=-200; dt=1; vx=.4416; vy=0; x1=0; y1=0;
	// radial distance to 3rd power
    r3=Math.pow((x-200)*(x-200)+(y-150)*(y-150), 1.5);
	// initial time step for velocity
    vx=vx+(x-200)/r3*dt/2; vy=vy+(y-150)/r3*dt/2;
    phase=1; 
    do {
	  // step on coordinates 
      x=x+vx*dt; y=y+vy*dt;   
      // radial distance to 3rd power	  
      r3=Math.pow((x-200)*(x-200)+(y-150)*(y-150), 1.5);
	  // step on velocities
      vx=vx+(x-200)/r3*dt; vy=vy+(y-150)/r3*dt;
	  // draw piece of curve
      if (phase==2) {  // continue drawing until out of scope
        x2=(int)(x+.5); y2=(int)(y+.5);
        g.drawLine(x1, y1, x2, y2); x1=x2; y1=y2;
        if ((x2>400)||(y2>300)) {
          phase=3;
        }
      }
      if ((phase==1)&&(x>=0)) {
        x1=(int)(x+.5); y1=(int)(y+.5); phase=2;
      }
    } while ((phase<3)&&(r3<66000000));
    return Math.atan2(vy, vx)/Math.PI*180;  // scattering angle
  }
  double ctr(double y) {        
	// Feynman Algorithm. No trajectory
    double x, vx, vy, dt, r3;   
    x=-200; dt=1; vx=.4416; vy=0;
    // radial distance to 3rd power	  
    r3=Math.pow((x-200)*(x-200)+(y-150)*(y-150), 1.5);
    // step on velocities
    vx=vx+(x-200)/r3*dt/2; vy=vy+(y-150)/r3*dt/2;
    do {
	  // step on coordinates and velocities 
      x=x+vx*dt; y=y+vy*dt;
      r3=Math.pow((x-200)*(x-200)+(y-150)*(y-150), 1.5);
      vx=vx+(x-200)/r3*dt; vy=vy+(y-150)/r3*dt;
    } while (r3<66000000);
    return Math.abs(Math.atan2(vy, vx)/Math.PI*180);  // scattering angle
  }
  double fcn(double t) {                               
	// Analytic N(theta)
    double z;
    z=Math.tan(Math.PI/180*t/2);
    return 10*Math.PI/9*(1+z*z)/z/z/z;
  }

  double rnd() {                                       
	// Random impact parameter
    double x, y, r;           
    int test=0;
    do {
      x=Math.random(); y=Math.random();
      // radial distance
      r=Math.pow((x-.5)*(x-.5)+(y-.5)*(y-.5), .5);
      if (r<=.5) { test=1; }
    } while (test==0);
    r=50*r;   // adjusting for screen size in pixels
    if (y<.5) { r=-r; }
    return 150+r;    
  }
  public void paint(Graphics g) {    // Paint
    double t; double z=125;
    int x1, y1, x2, y2, xorg, yorg;
    int [] cc=new int[20];
    if (kk==0) {  // button bounds
      b1.setBounds(360, 110, 30, 20);
      b2.setBounds(360, 140, 30, 20);
      b3.setBounds(360, 170, 30, 20);
    }
	// applet frame
    g.setColor(Color.black); g.drawRect(0, 0, 400, 300);
    if (kk<3) {    // draw scattering center
      g.setColor(Color.green);
      g.fillArc(190, 140, 20, 20, 0, 360);
      g.setColor(Color.blue);
      if (kk==1) {      // draw stepped impact parameter curves
        while (z<176) {                            
          hit(z, g);
          z=z+2;
        }
      }
      if (kk==2) {       // Random distribution
        for (int i=0; i<=25; i=i+1) {
          hit(rnd(), g);
        }
      }
    }
    if (kk==3) {       // 500 particle counter
      g.setFont(new Font("TimesRoman", Font.PLAIN, 14));
	  // draw frame
      xorg=65; yorg=260;
      g.setColor(Color.black); g.drawRect(xorg, 20, 270, 240);
      g.drawString("0", xorg-3, yorg+20);
      g.drawString("180", 335-9, yorg+20);
      g.drawString("Theta  deg", 160, yorg+20);
      g.drawString("Counts", 20, 146);
      g.setColor(Color.blue); g.drawString("0", xorg-12, yorg+6);
      g.drawString("480", xorg-24, 26);
      g.setColor(Color.magenta); g.drawString("0", 335+6, yorg+6);
      g.drawString("48", 335+6, 26);
	  // Analytic function
      g.setColor(Color.blue);      
      t=2*Math.atan(.2)/Math.PI*180;
      x1=xorg+(int)(1.5*t+.5); y1=yorg-(int)(fcn(t)/2+.5);
      g.drawLine(x1, yorg, x1, y1);
      for (int i=24; i<=60; i=i+2) {   // Range 480
        x2=xorg+3*i/2; y2=yorg-(int)(fcn((double)i)/2+.5);
        g.drawLine(x1, y1, x2, y2); x1=x2; y1=y2;
      }
      g.setColor(Color.magenta); y1=yorg-(int)(5*fcn(60.)+.5);
      for (int i=60; i<=180; i=i+2) {    // Range 48
        x2=xorg+3*i/2; y2=yorg-(int)(5*fcn((double)i)+.5);
        g.drawLine(x1, y1, x2, y2); x1=x2; y1=y2;
      }
	  // histograms
      g.setColor(Color.red);
      for (int i=1; i<=500; i=i+1) {
        x1=(int)(ctr(rnd())/10); cc[x1]=cc[x1]+1; x2=xorg+15*x1;
        if (x1<6) {
          y2=yorg-cc[x1]/2; g.drawLine(x2, y2, x2+15, y2);
        }
        else {
          y2=yorg-5*cc[x1]; g.fillRect(x2, y2, 15, 5);
        }
      }
      x1=0; x2=0;
      for (int i=2; i<=5; i=i+1) {
        x1=x1+cc[i];
      }
      for (int i=6; i<=17; i=i+1) {
        x2=x2+cc[i];
      }
	  // counting of particles
      g.setColor(Color.blue);
      g.drawString("N = "+x1, 265, 60);   // N with theta<60 deg
      g.setColor(Color.magenta);
      g.drawString("N = "+x2, 265, 80);   // N with theta>60 deg 
    }
    kk=0;
  }

  public void actionPerformed(ActionEvent e) {   // Action Buttons
    String tst;
    tst=e.getActionCommand();
    if (b1s.equals(tst)) { kk=1; }  // stteped impact parameter
    if (b2s.equals(tst)) { kk=2; }  // random impact parameter
    if (b3s.equals(tst)) { kk=3; }  // counts
    repaint();
  }
}