#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <GL/glut.h>
#include <GL/glu.h>
#include <GL/gl.h>
#include <assert.h>
#include <math.h>

// during initialization, initialize glut with GL_DOUBLE_BUFFER, as before.

#define IMAGE_MODE  GL_ABGR_EXT      // fastest pixel ordering
#define RES 500

int windowWidth = RES, windowHeight = RES;
int pass = 0;

void render();
void setColor(int r, int g, int b);

GLubyte *buffer;

void Redraw(void)
{
   render();
}

void Idle(void)
{
   glutPostRedisplay();
}

void rerender() {
   pass = 0;
   glutPostRedisplay();
}

double xOrigin = 0;
double yOrigin = 0;
double scale = .25;

int N = 0;
int *sx, *sy, *ss;
void initSxy(int _N) {
   N = _N;
   sx = new int[N*N];
   sy = new int[N*N];
   ss = new int[N*N];
   int b;
   for (int n = N*N-1 ; n >= 0 ; n--) {
      sx[n] = sy[n] = ss[n] = 0;
      for (int i = N/2 ; i > 0 ; i /= 2) {
         b = n/i/i+1;
         sx[n] += ((b/2) % 2) * N/i/2;
         sy[n] += ((b-1) % 2) * N/i/2;
         ss[n] = b==1 ? N/i : ss[n] > 1 ? ss[n] : 1;
      }
   }
}

void initBuffer() {
    if (buffer)
       free(buffer);
    buffer = new GLubyte[windowWidth*windowHeight*4];
    for (N = 1 ; N * 64 < windowWidth ; N *= 2)
       ;
    initSxy(N);
    pass = 0;
}

void Reshape(int width, int height)
{
    windowWidth = width;
    windowHeight = height;
    initBuffer();
    glViewport(0, 0, windowWidth, windowHeight);
}

void Button(int button, int state, int x, int y)
{
   double X = (double)x / windowWidth;
   double Y = (double)y / windowWidth;
   switch (button) {
   case GLUT_LEFT_BUTTON:
      if (state == GLUT_DOWN) {
         scale *= 2;
         xOrigin += xOrigin - X;
         yOrigin += yOrigin - Y;
         rerender();
      }
      break;
   case GLUT_MIDDLE_BUTTON:
      if (state == GLUT_DOWN)
        printf("%g %g %g\n", xOrigin, yOrigin, scale);
      break;
   case GLUT_RIGHT_BUTTON:
      if (state == GLUT_DOWN) {
         xOrigin += (X - xOrigin) / 2;
         yOrigin += (Y - yOrigin) / 2;
         scale /= 2;
         rerender();
      }
      break;
   }
}

int limit = 256;

void Key(unsigned char key, int x, int y)
{
    switch (key) {
    case 27:
       exit(0);
    case ' ':
       xOrigin = 0;
       yOrigin = 0;
       scale = .25;
       rerender();
       break;
    }
    if (key >= '0' && key <= '9') {
       limit = (32 << (key - '0'));
       rerender();
    }

}

void main(int argc, char *argv[])
{
   glutInit(&argc, argv);
   glutInitDisplayMode(GLUT_RGB);
   glutInitWindowSize(windowWidth, windowHeight);
   glutCreateWindow("Mandelbrot");
   
   // set up world space to screen mapping
   glMatrixMode(GL_PROJECTION);
   glLoadIdentity();
   glOrtho(0.0f, (GLfloat) windowWidth, (GLfloat) windowHeight, 0.0f, -1.0f, 1.0f);
   glMatrixMode(GL_MODELVIEW);
   glLoadIdentity();
   glViewport(0, 0, windowWidth, windowHeight);
   
   glClearColor(1.0,.95,.85,0);
   glutDisplayFunc(Redraw);
   glutIdleFunc(Idle);
   glutKeyboardFunc(Key);
   glutMouseFunc(Button);
   glutReshapeFunc(Reshape);

   glutMainLoop();
}

void setColor(int r, int g, int b) { glColor3ub(r, g, b); }
void fillRect(int x, int y, int w, int h) {
    glBegin(GL_QUADS);
       glVertex2i(x, y);
       glVertex2i(x, y+h);
       glVertex2i(x+w, y+h);
       glVertex2i(x+w, y);
    glEnd();
}

double color[3];

double mandelbrot(double x, double y) {
   double i,j,k,r;
   for(i=j=k=r=0;j=r*r-i*i+x,i=2*r*i+y,i*i+j*j<11&&k++<limit; r=j);
   return (double)k/limit;
}

void computeColor(double x, double y, double color[3]) {
   double t = mandelbrot(
      (x - xOrigin)/scale - 2.5,
      (y - yOrigin)/scale - 2.*windowHeight/windowWidth);
   t = (t > 1) ? 0 : t / 4;
   if (t < 1./16)
      t *= 16;
   else {
      t = 1 - (t-1./16) * 16./15;
      t*=t; t*=t; t*=t; t*=t; t*=t;
   }
   color[0] = t;
   t = t * t;
   color[1] = .85*t;
   t = t * t;
   color[2] = .6*t;
}

void render() {
   if (pass >= N*N)
      return;

   glDrawBuffer(GL_BACK); // draw into the back buffer

   int s = ss[pass];
   for (int y = sy[pass] ; y < windowHeight ; y += N)
   for (int x = sx[pass] ; x < windowWidth  ; x += N) {
      computeColor((double)x/RES, (double)y/RES*windowHeight/windowWidth, color);
      setColor((int)(255*color[0]), (int)(255*color[1]), (int)(255*color[2]));
      fillRect(x-s/2,y-s/2,s,s);
   }
   pass++;

   // copy the back buffer to the front:

   glReadPixels(0,0,windowWidth,windowHeight,IMAGE_MODE,GL_UNSIGNED_BYTE, buffer);
   glDrawBuffer(GL_FRONT);
   glDrawPixels(windowWidth,windowHeight,IMAGE_MODE,GL_UNSIGNED_BYTE,buffer);
}