/*
 *
 *  Este programa implementa a rasterizacao de 
 *  segmentos e poligonos por, respectivamente,
 *  algoritmo de Bresenham (tecla b) e de scanline 
 *  (tecla s). Além disso,
 *  ilustra a transformada discreta de Fourier 
 *  da imagem rasterizada e sua inversa. O programa
 *  é iniciado com uma imagem rasterizada de um
 *  segmento horizontal.
 *
 *            Ting (20/07/08)
 */

//
#define _USE_MATH_DEFINES

#include <stdio.h>      // Header file for standard file i/o.
#include <stdlib.h>     // Header file for malloc/free.
#include <math.h>       // Header file for math functions.
#include <GL/glut.h>    // Header File For The GLUT Library

static GLint height;

/* ascii code for the escape key */
#define ESCAPE 27

/* dimensoes da imagem */
#define ALTURA 64
#define LARGURA 64
#define scale 3
#define EPSILON 1.E-1
#define TRUE 1
#define FALSE 0

typedef struct Point2D_struct {   /* pontos */
  GLfloat x, y;
} Point2D;

typedef struct Y_bucket_struct {  /* y-bucket */
        GLint y_max;
        GLfloat x, m;
        struct Y_bucket_struct *proximo;
} Y_bucket;

unsigned int c = 1; //compressao da escala do espectro
int win;  // id da janela principal
int swin1, swin2, swin3; // id das subjanelas do glut
float *fourier, *Re, *Im;
GLubyte *ift_imagem, *ft_imagem, *s_imagem, *s_ft_imagem, *s_ift_imagem;
GLubyte imagem[LARGURA][ALTURA];

Y_bucket *TabelaDeArestasAtivas[ALTURA];  /* tabela de arestas ativas */

/************************************************************************/
/* Transformada Discreta de Fourier */

void dft(void) {
  int u, v, x, y;
  float t, t1, aux;

  float max=-1.E12, min=1E12;

  printf("Iniciou Transformada Discreta de Fourier\n");

  t = 2*M_PI;
  t1 = 1/((float)LARGURA * (float)ALTURA);
  for (u=0; u<LARGURA; u++) {
    for (v=0; v<ALTURA; v++) {
      Re[u*ALTURA+v] = Im[u*ALTURA+v] = 0.0;
      for (x=0; x<LARGURA; x++) 
	for (y=0; y<ALTURA; y++) {
	  aux = -t * (((float)u*x)/(float)ALTURA+((float)v*y)/(float)LARGURA);
            Re[u*ALTURA+v] += t1*(float)imagem[x][y] * cos (aux);
            Im[u*ALTURA+v] += t1*(float)imagem[x][y] * sin (aux);

	}
      // Obter a magnitude da transformada 
      fourier[u*ALTURA+v] = sqrt(Re[u*ALTURA+v]*Re[u*ALTURA+v] + Im[u*ALTURA+v]*Im[u*ALTURA+v]);
      if (fourier[u*ALTURA+v] > max) max = fourier[u*ALTURA+v];
      if (fourier[u*ALTURA+v] < min) min = fourier[u*ALTURA+v];
    }
  }

  /* Normalizar em relacao ao valor maximo (c=0) ou comprimir os
     valores logaritmicamente (c=1); e deslocar os valores para
     exibi-los de forma simetrica
  */
  if (c == 0) {
    for (u=0; u<LARGURA; u++) {
      for (v=0; v<ALTURA; v++) {
	ft_imagem[((u+ALTURA/2)%ALTURA)*ALTURA+(v+LARGURA/2)%LARGURA] = (GLubyte)(255.*fourier[u*ALTURA+v]/max);
      }
    }
  } else {
    max = log(1+max);
    for (u=0; u<LARGURA; u++) {
      for (v=0; v<ALTURA; v++) {
	ft_imagem[((u+ALTURA/2)%ALTURA)*ALTURA+(v+LARGURA/2)%LARGURA] = (GLubyte)(255.*log(1+fourier[u*ALTURA+v])/max);
      }
    }
  }
  printf ("Terminou TDF\n");
}

/* Transformada Discreta de Fourier inversa */

void idft(void) {
  int u, v, x, y;
  float t, aux, real, img;

  printf("Iniciou Transformada Discreta de Fourier Inversa\n");

  t = 2*M_PI;
  for (y=0; y<ALTURA; y++) {
    for (x=0; x<LARGURA; x++) {
      real = img = 0.0;
      for (v=0; v<ALTURA; v++)	
        for (u=0; u<LARGURA; u++) {
	  aux = t * (((float)v*y)/(float)ALTURA+((float)u*x)/(float)LARGURA);
          real += (float)(Re[u*ALTURA+v]*cos (aux)-Im[u*ALTURA+v]*sin(aux));
          img += (float)(Re[u*ALTURA+v]*sin (aux)+Im[u*ALTURA+v]*cos(aux));

	}
      if (fabs(img) > EPSILON) {
	printf("Erro na inversao: img = %lf\n", img);
      } 
      ift_imagem[x*ALTURA+y] = sqrt(real*real + img*img);
    }
  }

  printf("Terminou TDFI\n");

}

/************************************************************************/
// Rasterizar
// Linhas: y; Colunas:x ([y][x])

void makeBresImage(int x0, int y0, int x1, int y1)
{
   int y, x, inc1, inc2;
   float d;

   for (y = 0; y < ALTURA; y++) {
      for (x = 0; x < LARGURA; x++) {
         imagem[y][x] = (GLubyte) 0;
      }
   }

   if (y0 == y1) {
     for (x=x0; x< (x1+1); x++) 
      imagem[y0][x] = (GLubyte) 255;
     return;
   }

   if (x0 == x1) {
     for (y=y0; y< (y1+1); y++) 
      imagem[y][x0] = (GLubyte) 255;
     return;
   }

   y = y0;
   if (x1 > x0 && y1 > y0 && (x1-x0) >= (y1-y0)) {
     d = 2*(y0-y1)*(x0+1)+(x1-x0)*(2*y0+1)+2*x0*y1 - 2*x1*y0;
     inc1 = 2*(x1-x0)+2*(y0-y1);
     inc2 = 2*(y0-y1);

     for (x=x0; x< (x1+1); x++) {
      imagem[y][x] = (GLubyte) 255;
      if (d < 0) {
	y = y+1;
	d = d+inc1;
      } else {
	d = d+inc2;
      }
     }
   } else if (x1 > x0 && y1 > y0 && (y1-y0) > (x1-x0)){
     d = 2*(x0-x1)*(y0+1)+(y1-y0)*(2*x0+1)+2*y0*x1 - 2*y1*x0;
     inc1 = 2*(y1-y0)+2*(x0-x1);
     inc2 = 2*(x0-x1);

     for (y=y0; y < (y1+1); y++) {
      imagem[y][x] = (GLubyte) 255;
      if (d < 0) {
	x = x+1;
	d = d+inc1;
      } else {
	d = d+inc2;
      }
     }
   } else if (x0 > x1 && y1 > y0 && (y1-y0) >= (x0-x1)) {
     d = 2*(x1-x0)*(y0+1)+(y1-y0)*(2*x1+1)+2*y0*x0 - 2*y1*x1;
     inc1 = 2*(y1-y0)+2*(x1-x0);
     inc2 = 2*(x1-x0);

     for (y=y0; y < (y1+1); y++) {
      imagem[y][x] = (GLubyte) 255;
      if (d < 0) {
	x = x-1;
	d = d+inc1;
      } else {
	d = d+inc2;
      }
     }
   } 
}

// Scanline: construcao de y-bucket
void inserirYbucket(GLint i, Y_bucket *y_bucket)
{
  Y_bucket *pre, *tmp;

  /* Determinar y_bucket que contenha x maior que y_bucket->x */
  for (pre=0, tmp=TabelaDeArestasAtivas[i];
           tmp && tmp->x < y_bucket->x;
           pre=tmp, tmp = tmp->proximo);
  /* Inserir y_bucket na tabela de arestas ativas */
  if (pre) {
        y_bucket->proximo = tmp;
        pre->proximo = y_bucket;
  } else {
        y_bucket->proximo = TabelaDeArestasAtivas[i];
        TabelaDeArestasAtivas[i] = y_bucket;
  }
}

void makeScanImage(int n, Point2D *vertices)
{
  GLint i, j, p;
  GLfloat y_max, x, delta_x, delta_y, x_1, x_2, m;
  Y_bucket *y_bucket, *tmp;

  // Fundo Preto
  for (j = 0; j < ALTURA; j++) {
      for (i = 0; i < LARGURA; i++) {
         imagem[i][j] = (GLubyte) 0;
      }
  }

  /* Inicializar a tabela de arestas ativas */
  for (i=0; i<n; i++) {
    /* Ordenar a orientacao das arestas de forma que y seja sempre
       crescente */
    if (i==(n-1)) p = 0;
    else p =i+1;

    if (vertices[p].y > vertices[i].y) {
      j = vertices[i].y;
      y_max = vertices[p].y;
      x = vertices[i].x;
      delta_y = vertices[p].y - vertices[i].y;
      delta_x = vertices[p].x - vertices[i].x;
      m = delta_x/delta_y;
    } else {
      j = vertices[p].y;
      y_max = vertices[i].y;
      x = vertices[p].x;
      delta_y = vertices[i].y - vertices[p].y;
      delta_x = vertices[i].x - vertices[p].x;
      m = delta_x/delta_y;
    }
    /* Construir y_bucket */
    y_bucket = (Y_bucket *)malloc(sizeof(Y_bucket));
    y_bucket->y_max = y_max;
    y_bucket->x = x;
    y_bucket->m = m;
    y_bucket->proximo = 0;

    /* Inserir na lina j da tabela de arestas ativas */
    if (TabelaDeArestasAtivas[j] == 0)
      TabelaDeArestasAtivas[j] = y_bucket;
    else {
      /* Determinar y_bucket que contenha x maior que y_bucket->x */
      inserirYbucket(j, y_bucket);
    }
  }

  /* Rasterizar */
  for (i=0; i < ALTURA; i++) {
    if (TabelaDeArestasAtivas[i]) {
      do {
        tmp = TabelaDeArestasAtivas[i];
        x_1 = tmp->x;
        x_2 = tmp->proximo->x;
        for (j=x_1; j < x_2; j++)
	  imagem[i][j] = (GLubyte) 255;;

        /* Avancar para o proximo intervalo da linha */
        TabelaDeArestasAtivas[i] = tmp->proximo->proximo;
        if ((i+1) < tmp->proximo->y_max) {
          /* Atualizar o valor de x */
          tmp->proximo->x = tmp->proximo->x + tmp->proximo->m;
          /* Transferir y_buckets para a proxima linha de varredura */
          inserirYbucket(i+1, tmp->proximo);
        } else {
          /* Remover y_buckets */
          free(tmp->proximo);
        }
        if ((i+1) < tmp->y_max) {
          /* Atualizar o valor de x */
          tmp->x = tmp->x + tmp->m;
          /* Transferir y_buckets para a proxima linha de varredura */
          inserirYbucket(i+1, tmp);
        } else {
          /* Remover y_buckets */
          free(tmp);
        }
      } while (TabelaDeArestasAtivas[i]);
    }
  }
}

/************************************************************************/
// Miscellaneous functions

void update(void)
{
  dft();
  idft();
  glutSetWindow(swin1);
  glutPostRedisplay();
  glutSetWindow(swin2);
  glutPostRedisplay();
  glutSetWindow(swin3);
  glutPostRedisplay();
}

/************************************************************************/
// OpenGL functions
//
void init(void)
{
     
  /* Fourier Transform */
  fourier = (float *) malloc(sizeof(float)*LARGURA*ALTURA);
  if (fourier == NULL) {
	printf("Error allocating memory for fourier transform data");
        exit(0);
	return;
  }
  Re = (float *) malloc(sizeof(float)*LARGURA*ALTURA);
  if (Re == NULL) {
	printf("Error allocating memory for fourier transform data");
        exit(0);
	return;
  }
  Im = (float *) malloc(sizeof(float)*LARGURA*ALTURA);
  if (Im == NULL) {
	printf("Error allocating memory for fourier transform data");
        exit(0);
	return;
  }
  ft_imagem = (GLubyte *) malloc(sizeof(GLubyte)*LARGURA*ALTURA);
  if (ft_imagem == NULL) {
	printf("Error allocating memory for fourier transform image data");
        exit(0);
	return;
  }

  ift_imagem = (GLubyte *) malloc(sizeof(GLubyte)*LARGURA*ALTURA);
  if (ift_imagem == NULL) {
	printf("Error allocating memory for fourier transform image data");
        exit(0);
	return;
  }

  /* Space for scaled image */
  s_imagem = (GLubyte *) malloc(sizeof(GLubyte)*LARGURA*ALTURA*scale*scale);
  if (s_imagem == NULL) {
	printf("Error allocating memory for scaled imagem data");
        exit(0);
	return;
  }

  s_ft_imagem = (GLubyte *) malloc(sizeof(GLubyte)*LARGURA*ALTURA*scale*scale);
  if (s_ft_imagem == NULL) {
	printf("Error allocating memory for fourier transform image data");
        exit(0);
	return;
  }

  s_ift_imagem = (GLubyte *) malloc(sizeof(GLubyte)*LARGURA*ALTURA*scale*scale);
  if (s_ft_imagem == NULL) {
	printf("Error allocating memory for fourier transform image data");
        exit(0);
	return;
  }

   glClearColor (0.0, 0.0, 0.0, 0.0);
   glShadeModel(GL_FLAT);
   glPixelStorei(GL_UNPACK_ALIGNMENT, 1);

   /* Rasterize um segmento horizontal */
   makeBresImage(4,16,60,16);
   /* transformada discreta */
   dft();
   /* transformada inversa */
   idft();
}

void display_swin1()
{
  // Desenhar os pixels da imagem
   glutSetWindow(swin1);
   glDisable(GL_DITHER);

   glClear(GL_COLOR_BUFFER_BIT);
   glRasterPos2i(0, 0);
   gluScaleImage(GL_LUMINANCE, LARGURA, ALTURA, GL_UNSIGNED_BYTE,
                 imagem,scale*LARGURA, scale*ALTURA, GL_UNSIGNED_BYTE,s_imagem);
   glDrawPixels(scale*LARGURA, scale*ALTURA, GL_LUMINANCE, GL_UNSIGNED_BYTE, s_imagem);
   glFlush();
}

void display_swin2()
{
  // Desenhar os pixels da imagem da transformada
  glutSetWindow(swin2);
  glDisable(GL_DITHER);
  glClear(GL_COLOR_BUFFER_BIT);
  glRasterPos2i(0, 0);
  gluScaleImage(GL_LUMINANCE, LARGURA, ALTURA, GL_UNSIGNED_BYTE,
                 ft_imagem,scale*LARGURA, scale*ALTURA, GL_UNSIGNED_BYTE,s_ft_imagem);
  glDrawPixels(scale*LARGURA, scale*ALTURA, GL_LUMINANCE, GL_UNSIGNED_BYTE, s_ft_imagem);
  glFlush();
}

void display_swin3()
{
  // Desenhar os pixels da imagem da transformada
  glutSetWindow(swin3);
  glDisable(GL_DITHER);
  glClear(GL_COLOR_BUFFER_BIT);
  glRasterPos2i(0, 0);
  glDrawPixels(ALTURA,LARGURA,GL_LUMINANCE, GL_UNSIGNED_BYTE, ift_imagem);

  gluScaleImage(GL_LUMINANCE, LARGURA, ALTURA, GL_UNSIGNED_BYTE,
                 ift_imagem,scale*LARGURA, scale*ALTURA, GL_UNSIGNED_BYTE,s_ift_imagem);
  glDrawPixels(scale*LARGURA, scale*ALTURA, GL_LUMINANCE, GL_UNSIGNED_BYTE, s_ift_imagem);
  glFlush();
}

void
main_display(void)
{
    glClearColor(0.8, 0.8, 0.8, 0.0);
    glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
    glColor3f(0.f, 0.f, 0.f);
	display_swin1();
	display_swin2();
	display_swin3();
}

void main_reshape(int w, int h)
{
   glViewport(0, 0, (GLsizei) w, (GLsizei) h);
   height = (GLint) h;
   glMatrixMode(GL_PROJECTION);
   glLoadIdentity();
   gluOrtho2D(0.0, (GLdouble) w, 0.0, (GLdouble) h);
   glMatrixMode(GL_MODELVIEW);
   glLoadIdentity();

  /* redefinir a posição e as dimensões da sub-janela 1 */
  glutSetWindow(swin1);
  glutPositionWindow(0, 0);
  glutReshapeWindow(scale*LARGURA, scale*ALTURA);
  /* redefinir a posição e as dimensões da sub-janela 2 */
  glutSetWindow(swin2);
  glutPositionWindow(scale*LARGURA, 0);
  glutReshapeWindow(scale*LARGURA, scale*ALTURA);
  /* redefinir a posição e as dimensões da sub-janela 3 */
  glutSetWindow(swin3);
  glutPositionWindow(2*scale*LARGURA, 0);
  glutReshapeWindow(scale*LARGURA, scale*ALTURA);
}

void reshape_swin1(int w, int h)
{
  glMatrixMode(GL_PROJECTION);
  glLoadIdentity();
  gluOrtho2D(0.0, (GLfloat) scale*LARGURA, 0.0, (GLfloat) scale*ALTURA);
  glMatrixMode(GL_MODELVIEW);
  glLoadIdentity();
}

void reshape_swin2(int w, int h)
{
  glMatrixMode(GL_PROJECTION);
  glLoadIdentity();
  gluOrtho2D(0.0, (GLfloat) scale*LARGURA, 0.0, (GLfloat) scale*ALTURA);
  glMatrixMode(GL_MODELVIEW);
  glLoadIdentity();
}

void reshape_swin3(int w, int h)
{
  glMatrixMode(GL_PROJECTION);
  glLoadIdentity();
  gluOrtho2D(0.0, (GLfloat) scale*LARGURA, 0.0, (GLfloat) scale*ALTURA);
  glMatrixMode(GL_MODELVIEW);
  glLoadIdentity();
}

void keyboard(unsigned char key, int x, int y)
{
  //   Point2D pontos[4] = {
  //   {16.0, 13.0}, {13.0, 16.0}, {45.0, 48.}, {48.0, 45.}};
   Point2D pontos[4] = {
     {16.0, 13.0}, {55.0, 16.0}, {55.0, 48.}, {10.0, 32.}};

   switch (key) {
   case 'c':
     if (c==0) c=1;
     else c=0;
     dft();
     glutSetWindow(swin2);
     glutPostRedisplay();
     break;
   case 'b':
     makeBresImage(4,16,60,60);
     update();
     break;
   case 's':
     makeScanImage(4,pontos);
     update();
     break;
   case 27:
     free((char *)ft_imagem);
     free((char *)ift_imagem);
     free((char *)fourier);
     free((char *)Re);
     free((char *)Im);
     free((char *)s_ft_imagem);
     free((char *)s_ift_imagem);
     free((char *)s_imagem);
         exit(0);
         break;
      default:
         break;
   }
}

/************************************************************************/
// Main configuration
//
int main(int argc, char** argv)
{

   glutInit(&argc, argv);
   glutInitDisplayMode(GLUT_SINGLE | GLUT_RGB );
   glutInitWindowSize(LARGURA*scale*3, ALTURA*scale);
   glutInitWindowPosition(0,0);
   win = glutCreateWindow(argv[0]);
    glutDisplayFunc(main_display);
    glutReshapeFunc(main_reshape);
    glutKeyboardFunc(keyboard);
   swin1 = glutCreateSubWindow(win,0,0,scale*LARGURA,scale*ALTURA);
    glutReshapeFunc(reshape_swin1);
    glutDisplayFunc(display_swin1);
    glutKeyboardFunc(keyboard);
   swin2 = glutCreateSubWindow(win,scale*LARGURA,0,scale*LARGURA,scale*ALTURA);
    glutReshapeFunc(reshape_swin2);
    glutDisplayFunc(display_swin2);
    glutKeyboardFunc(keyboard);
   swin3 = glutCreateSubWindow(win,2*scale*LARGURA,0,scale*LARGURA,scale*ALTURA);
    glutReshapeFunc(reshape_swin3);
    glutDisplayFunc(display_swin3);
    glutKeyboardFunc(keyboard);

   init();

   glutMainLoop();
   return 0; 
}