#include <iostream.h>
#include <math.h>
#include <GL/gl.h>
#include <GL/glu.h>
#include <GL/glut.h>

//-------------------------------------------------------------------------------
//Polar to Rectangular Coodinate Conversion Functions (and Inverses)

double PtoRx(double r,double t)
  { return r*cos(t); }

double PtoRy(double r,double t)
  { return r*sin(t); }

double RtoPr(double x,double y)
  { return sqrt(x*x+y*y); }

double RtoPt(double x,double y)
  { return atan(y/x); }

//------------------------------------------------------------------------------
// OpenGL Calls

const XWINDOW=500;
const YWINDOW=500;

const C_XMin=0;
const C_XMax=500;
const C_YMin=0;
const C_YMax=500;

const TargetX=100;
const TargetY=0;

double r=100,t=30*M_PI/180;

static void Init()
  { glClearColor(0.0,0.0,0.0,1.0);
    glShadeModel(GL_FLAT);
  }

static void Reshape(int w,int h)
  { glViewport(0,0,(GLint)w,(GLint)h);
    glMatrixMode(GL_PROJECTION);
    glLoadIdentity();
    glOrtho(C_XMin,C_XMax,C_YMin,C_YMax,-1.0,1.0);
    glMatrixMode(GL_MODELVIEW);
  }

int InBounds(double x,double y)
  { return ((x>=C_XMin)&&(x<=C_XMax)); }

int AlignedHoriz(double x,double y)
  { return (y==TargetY); }

int AlignedVert(double x,double y)
  { return (x==TargetX); }

double DistanceFromTarget(double x,double y)
  { return sqrt((TargetX-x)*(TargetX-x)+(TargetY-y)*(TargetY-y)); }

double simulate(double a,double b)
  { double x=0,y=0,q=0,d=DistanceFromTarget(0,0);
    while(y>=TargetY)
        { d=DistanceFromTarget(x,y);
          x=q*a;
          y=q*b-16*q*q;
          q+=0.1;
        }
    return d;
  }

void Display()
  { glFlush();
    glClear(GL_COLOR_BUFFER_BIT);

    double x,y,t;
    glBegin(GL_POINTS);
      for(y=0;y<YWINDOW;y+=1)
        { for(x=0;x<XWINDOW;x+=1)
            { t=simulate(C_XMin+x*(C_XMax-C_XMin)/XWINDOW,
                         PtoRy(C_XMin+x*(C_XMax-C_XMin)/XWINDOW,y/YWINDOW*M_PI));
              glColor3f(1-sqrt(t)/10,0.0,0.0);
              
              //cout<<"x= "<<x<<" y= "<<y<<" sqrt(t)/1000= "<<sqrt(t)/10<<"\n";
              glVertex3f(C_XMin+x*(C_XMax-C_XMin)/XWINDOW+0.5,C_YMin+y*(C_YMax-C_YMin)/YWINDOW+0.5,0);
            } 
        }
   glEnd();


  }

void KeyboardInput(unsigned char key,int x,int y)
  { switch (key)
      { case 27: //ESC Quits
          exit(0);
          break;
      }
  }

void MouseInput(int button,int state,int x,int y)
  { switch (button)  
      { case (GLUT_LEFT_BUTTON):    
          { if (state==GLUT_DOWN)
              { cout<<"x= "<<C_XMin+x*(C_XMax-C_XMin)/XWINDOW
                    <<" y= "<<C_YMin+(YWINDOW-y)*(C_YMax-C_YMin)/YWINDOW
                    <<" sim= "
	            <<simulate(C_XMin+x*(C_XMax-C_XMin)/XWINDOW,
                               C_YMin+(YWINDOW-y)*(C_YMax-C_YMin)/YWINDOW)
                    <<"\n"; 
              }
          }    
          break;
      }
  }

int main(int argc,char **argv)
  { glutInit(&argc,argv);
    glutInitDisplayMode(GLUT_RGBA);
    glutInitWindowSize(XWINDOW,YWINDOW);
    glutInitWindowPosition(0,0);
    glutCreateWindow(argv[0]);
    Init();
    glutReshapeFunc(Reshape);
    glutDisplayFunc(Display);
    glutKeyboardFunc(KeyboardInput);
    glutMouseFunc(MouseInput); 
    glutMainLoop();
    return 0; 
  }