/***************************************************
* Neutron simulator program written by Vernon Brown
* in Quick Basic in 1991, version 911110 (date code)
* meaning November 10th 1991. This version copied
* by Vernon Brown in February 2006 using C++ and the 
* Linux SDL graphics library. The reason for writing
* the program was to determine how many wave cycles
* of a photon would occupy a resonant pattern to
* become one of the shells comprising a neutron.
* It was during this experiment in 1991 that Vernon
* Brown discovered the square-of-the-shells rule and
* the fundamental cause of the strong nuclear binding
* force.
* Requires SDL library and my personal library, veblib.h
* Compile with the console command: 
* g++ neut.cpp -l SDL
***************************************************/
#include <stdlib.h>
#include <SDL/SDL.h>
#include <stdio.h>
#include <math.h>
#include "veblib.h"

int hash360(int, double, unsigned int, int);
char version[] = {"20080705"};

/* global shell diameters to scale */
unsigned int gelec = 3497;
unsigned int gs1 = 700, gs2 = 274, gs3 = 42, gs4 = 1;

/*screenup screenleft screenright screendown  */
int sup = 0, slt=0, srt=800, sdn=800, wait=50;
int main(int argc, char *argv[]) {
    int dummy, x, y; 
    unsigned int color = 0xeeeeee;
    double rotate;
    double rotate2;
    double rotate3;
    double rotate4;
    int positive = 1, negative = 2;
    int ppo = 1;                    /* photon cycles per orbit    */
    char cmd[10] = {'1'};
    char c;
    if((SDL_Init(SDL_INIT_VIDEO|SDL_INIT_AUDIO)==-1)) { 
        printf("Could not initialize SDL: %s.\n", SDL_GetError());
        exit(-1);
    }
    dummy = init_screen(srt, sdn);
    int waiting = 10000;
    ppo = 1; x = 0; c = 0;
    dummy = ColorFill(color, srt, sdn);
    rotate = .1;
    while(waiting--)
    { 
       c = keystroke();
       ppo = atoi(cmd);
       rotate2 = rotate * 2;    /* Speed of rotation not to scale */
       rotate3 = rotate2 * 2 * PI;
       rotate4 = rotate3 * 2 * PI;
       /* Plot to scale in units of Shell Four diameter  */
       dummy = hash360(ppo, rotate,  gs1 / 2, negative); 
       dummy = hash360(ppo, rotate2, gs2 / 2, positive);
       dummy = hash360(ppo, rotate3, gs3 / 2, negative);
       dummy = hash360(ppo, rotate4, gs4, positive);  
       SDL_UpdateRect(screen, sup, slt, srt, sdn);
       SDL_Delay(wait);
       rotate += .1;
       if(wait > 10) wait--;
       dummy = ColorFill(color, srt, sdn);
       if(c == 'q') break;
       if(c == 'w') wait = 1000;
    }
    SDL_Quit();
    exit(0);
}


int hash360(int ppo, double rotate, unsigned int radius, int polarity)
  {
     int dummy;
     double rs = 2.549920405;  /* neutron - proton     */
     int mx = 20, ln = 0;      /* Hash Mark density    */
     int o1 = ppo;             /* Cycles per pattern   */

     unsigned int c1 = 0xff0000, c2 = 0x0000ff;
     unsigned int x, y, color;
     double aa = .001;   
     double a1 = aa * o1;  
     double ba = .001;
     double b1 = ba * o1;
     double mp, mn;
     double x1, x2, y1, y2;
     double j1;
     double k1;
     double maxamp = 40;
     double vall = maxamp;
     double v1 = vall, v2 = vall, v3 = vall, v4 = vall; /* Amplitude */
     double w, z, j, k, l, m;
     double rot;
     char c;
     rot = rotate / 6000;
     x1 = 0;

     if(polarity == 2){        /* 2 is negative 1 is positive */
       c1 = 0x0000ff;          /* Negative                    */
       c2 = 0xff0000;          /* Positive                    */
     }
     else{
       c1 = 0xff0000;
       c2 = 0x0000ff;
     }

     /* Amplitudes and spacing of hash marks */
     if(radius > 200)
       { maxamp = 50; mx = 10; }
     if(radius < 201 && radius > 100) 
       { maxamp = 50; mx = 5; }
     if(radius < 101 && radius > 30) 
       { maxamp = 10; mx = 3; }
     if(radius < 31)
       { maxamp = 5;  mx = 1;  }
     /***************************************/

     int makingcircle = 6280;
     while(makingcircle--)
        {    
           vall = abs(sin(x1 + rotate) * maxamp);
           v3 = vall, v4 = vall; /* Amplitude */

           x1 += aa;  /* increment degrees of rotation */
           y1 += ba;
           j1 += a1;
           k1 += b1;
            
            /* Calculate screen location to plot each dot */
           w = 400 + cos(x1) * radius;
           z = 400 + sin(y1) * radius;

           l = cos(j1) * v3;
           m = sin(k1) * v4;
            
            /* Plot each dot */

           x =  (unsigned int)w + (unsigned int)j;
           y =  (unsigned int)z - (unsigned int)k;
           color =  (unsigned int) c1;

            /* Make polarity hash marks  */

           if(ln < 1000) ln++;          /* space between hash marks */
           if(ln == mx){                /* determined by mx         */
              ln = 0;
              line(w, z, w+l, z+m, c1);
              line(w, z, w-l, z-m, c2);
           }
           if(x1 > 46) x1 = 1;
           if(y1 > 46) y1 = 1;
        }
        //SDL_UpdateRect(screen, sup, slt, srt, sdn);
     return 0;
  }