/* Pinball simulator */
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#define INF_T  10000.0
#define NO_HIT 100
#define ND 3
#define MAX_HIT 100
#define PI 3.14159265358979323846

/* Macros for doing vector arithmatic */
#define setv(v, x1, y1) {v.x = x1; v.y = y1;}
#define copyv(v, w) {v.x = w.x; v.y = w.y;}
#define dotv(v1, v2) (v1.x*v2.x + v1.y*v2.y)
#define hypotv(v) hypot(v.x, v.y)
#define smulv(v1, a, v2) {v1.x=a*v2.x; v1.y=a*v2.y;}
#define subv(v1, v2, v3) {v1.x=v2.x-v3.x; v1.y=v2.y-v3.y;}
#define addv(v1, v2, v3) {v1.x=v2.x+v3.x; v1.y=v2.y+v3.y;}
#define addvx(v1, v2, x1, y1) {v1.x=v2.x+x1; v1.y=v2.y+y1;}
#define addtov(v, x1, y1) {v.x=v.x+x1; v.y=v.y+y1;}
#define linear(v, a, va, b, vb) {v.x=a*va.x+b*vb.x; v.y=a*va.y+b*vb.y;}
#define polar(v, r, theta) {v.x=r*cos(theta); v.y=r*sin(theta);}
#define printv(v) printf("(%f, %f)\n", v.x, v.y)

typedef struct _Vector {
  double x;
  double y;
} Vector;

Vector disk[3];
double a=1.0/6.0, ia, asq;
char *progname;

main(argc, argv)		/* driver */
     int argc;
     char *argv[];
{
  int hits;
  double time, init_angle=0.0, l=0.5/sqrt(3.0);
  Vector pos, vel;
  int collisions();

  progname = argv[0];
  while (argc>2 && argv[1][0]=='-') {
    switch (argv[1][1]) {
    case 'a':
      a = atof(argv[2]);
      break;
    case 'R':
      a = 1.0/atof(argv[2]);
      break;
    default:
      fprintf(stderr, "usage: %s [-a x.xx|-R x.xx] angle\n", progname);
      exit(0);
      break;
    }
    argc-=2;
    argv+=2;
  }
  if (a>1.0/sqrt(3.0)) {
    fprintf(stderr, "error(%s): a=%f is too big\n", progname, a);
    exit(1);
  }
    
  ia = 1.0/a;
  asq = a*a;

  if (argc>1)
    init_angle = (PI/180.0)*atof(argv[1]);

  setv(disk[0], -0.5, -l);
  setv(disk[1], 0.5, -l);
  setv(disk[2], 0, 2.0*l);
  polar(vel, 1.0, init_angle);
  setv(pos, 0.0, 0.0);
  /* start at centre in the direction given by initial angle */
  hits = collisions(pos, vel, &time);
  printf("Number of hits = %d, time of flight = %f\n", hits, time);
}

/* calculate trajectory */
int collisions(pos, vel, total_time)
     Vector pos, vel;
     double *total_time;
{
  int len, i, cd, d_tmp, no_hits;
  Vector dd, perp;
  double tmp1, tmp2, t_tmp, time;

  printf("Initial position ");
  printv(pos);

  no_hits = 0;
  cd = NO_HIT;
  
  for (len=0; len<MAX_HIT; len++) {
    d_tmp = NO_HIT;
    time = INF_T;
    for(i=0; i<ND; i++) {
      if (i!=cd) {
	subv(dd, disk[i], pos);
	if ((tmp1 = dotv(vel, dd)) > 0.0) {
	  if ((tmp2=tmp1*tmp1-dotv(dd, dd)+asq) > 0.0) {
	    if ((t_tmp=tmp1-sqrt(tmp2))<time) {
	      d_tmp = i;
	      time = t_tmp;
	    }
	  }
	}
      }
    }
    if (time<INF_T) {
      no_hits++;
      *total_time += time;
      cd = d_tmp;
      linear(pos, 1.0, pos, time, vel);
      printf("collision[%d] with disk %d at time %f and position ", 
	     no_hits, cd+1, *total_time);
      printv(pos);
      linear(perp, ia, pos, -ia, disk[cd]);
      tmp1 = -2.0*dotv(vel, perp);
      linear(vel, 1.0, vel, tmp1, perp);
    } else {
      *total_time += time = 2.0;
      linear(pos, 1.0, pos, time, vel);
      printf("Final position ");
      printv(pos);
      break;
    }
  }
  return(no_hits);
}