path: root/Smoke/renderer_gl.c

#ifdef G_OS_WIN32
#define WIN32_LEAN_AND_MEAN 1
#include <windows.h>
#endif

#include <stdio.h>
#include <math.h>

#include <GL/gl.h>
#include <GL/glu.h>
#ifndef GL_VERSION_1_2
#include <gdk/glext/glext.h>
#endif

#include <rfftw.h>

#include "funcs.h"
#include "fluids.h"
#include "palette.h"
#include "normals.h"
#include "smoke.h"
#include "colormap.h"
#include "glyphs.h"
#include "isolines.h"
#include "streamlines.h"
#include "heightplots.h"
#include "flowvis.h"
#include "seedpoint.h"
#include "renderer_gl.h"

#define DEFAULT_X_POS 0.0f
#define DEFAULT_Y_POS 0.0f
#define DEFAULT_ZOOM -750.0f
#define MIN_ZOOM -100.0f
#define MAX_ZOOM -4000.0f
#define DEFAULT_ZOOM_SPEED 100

#define MIN_ZOOM -100.0f
#define MAX_ZOOM -4000.0f

#define LEGEND_X_POS -298.0f
#define LEGEND_Y_POS -295.0f
#define LEGEND_Z_POS -735.0f

#define BITMAP_FILESIZE 0x02
#define BITMAP_OFFSET 0x0a
#define BITMAP_HEADERSIZE 0x0e
#define BITMAP_WIDTH 0x12
#define BITMAP_HEIGHT 0x16

#define NORMAL_SCALE 16

float x_pos = DEFAULT_X_POS;
float y_pos = DEFAULT_Y_POS;
float z_pos = DEFAULT_ZOOM;

static int renderer_render_grid = FALSE;
static int renderer_zoomspeed = DEFAULT_ZOOM_SPEED;
static GLuint texture;
static int grid_cell_size = 50;



void renderer_load_texture(char *filename)
{
  FILE *bitmap;

  
  bitmap = fopen(filename, "rb");
  if (bitmap > 0) {
    unsigned int dataoffset, filesize;
    GLsizei width, height;
    unsigned char *imagedata;

    fseek(bitmap, BITMAP_FILESIZE, SEEK_SET);
    fread(&filesize, 4, 1, bitmap);
    fseek(bitmap, BITMAP_OFFSET, SEEK_SET);
    fread(&dataoffset, 4, 1, bitmap);
/*    fseek(bitmap, BITMAP_HEADERSIZE, SEEK_SET); */
/*    fread(&headersize, 4, 1, bitmap); */
    fseek(bitmap, BITMAP_WIDTH, SEEK_SET);
    fread(&width, 4, 1, bitmap);
/*    fseek(bitmap, BITMAP_HEIGHT, SEEK_SET); */
    fread(&height, 4, 1, bitmap);

    imagedata = (unsigned char *)malloc((size_t)(filesize -dataoffset) /* *(char *) */);

    fseek(bitmap, dataoffset, SEEK_SET);
    fread(imagedata, (size_t)(filesize -dataoffset), 1, bitmap);

    glGenTextures(1, &texture);
    glBindTexture(GL_TEXTURE_2D, texture);
    gluBuild2DMipmaps(GL_TEXTURE_2D, GL_RGB, width, height, GL_BGR, GL_UNSIGNED_BYTE, imagedata);
    glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MIN_FILTER, GL_LINEAR);
	  glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MAG_FILTER, GL_LINEAR);
    
    free(imagedata);
    fclose(bitmap);
  }
}


static void render_legend(void)
{
  int i, s, c;
  struct color4f color;

  glDisable(GL_LIGHTING);
  glDisable(GL_DEPTH_TEST);
  glDisable(GL_BLEND);

  glLoadIdentity();

  for (i = 0; i < winWidth; i++)
  {
		float value, clamp_scaled_min, clamp_scaled_max, scale_scaled_min, scale_scaled_max;

		value = (float)((float)i/winWidth);

    s = colormap_get_scaling();
    c = colormap_get_clamping();

    colormap_set_scaling(0);
    colormap_set_clamping(0);

    color = colormap_get_color(value);

    colormap_set_scaling(s);
    colormap_set_clamping(c);

		clamp_scaled_min = colormap_get_clamp_min() * winWidth;
		clamp_scaled_max = colormap_get_clamp_max() * winWidth;

		scale_scaled_min = colormap_get_scale_min() * winWidth;
		scale_scaled_max = colormap_get_scale_max() * winWidth;

		glBegin(GL_LINES);
			if (colormap_get_scaling()) {
				glColor3f(1, 0, 0);
				glVertex2i(scale_scaled_min, winHeight - 4);
				glVertex2i(scale_scaled_min, winHeight - 18);
				glVertex2i(scale_scaled_max, winHeight - 4);
				glVertex2i(scale_scaled_max, winHeight - 18);
			}
    	if (colormap_get_clamping()) {      
  			glColor4f(0.2, 0.2, 0.2, 1.0);
	  		glVertex2i(clamp_scaled_min, winHeight - 6);
		  	glVertex2i(clamp_scaled_min, winHeight - 18);
			  glVertex2i(clamp_scaled_max, winHeight - 6);
  			glVertex2i(clamp_scaled_max, winHeight - 18);
      }
			if (!(i %20)) {
				glColor3f(0.7, 0.7, 0.7);
				glVertex2i(i -1, winHeight - 6);
				glVertex2i(i -1, winHeight - 22);
				glVertex2i(i, winHeight - 6);
				glVertex2i(i, winHeight - 22);
				glVertex2i(i +1, winHeight - 6);
				glVertex2i(i +1, winHeight - 22);
			} else {
        glColor3f(color.r, color.g, color.b);
				glVertex2i(i, winHeight - 6);
				glVertex2i(i, winHeight - 18);
			}
		glEnd();

		glBegin(GL_TRIANGLES);
    	if (colormap_get_scaling()) {
				glColor3f(1, 0, 0);
				glVertex2i(scale_scaled_min, winHeight - 6);
				glVertex2i(scale_scaled_min - 4, winHeight);
				glVertex2i(scale_scaled_min + 4, winHeight);
				glVertex2i(scale_scaled_max, winHeight - 6);
				glVertex2i(scale_scaled_max - 4, winHeight);
				glVertex2i(scale_scaled_max + 4, winHeight);
			}
    	if (colormap_get_clamping()) {
			  glColor3f(0.9, 0.9, 0.9);
        glVertex2i(clamp_scaled_min, winHeight - 18);
        glVertex2i(clamp_scaled_min - 4, winHeight - 25);
        glVertex2i(clamp_scaled_min + 4, winHeight - 25);
        glVertex2i(clamp_scaled_max, winHeight - 18);
        glVertex2i(clamp_scaled_max - 4, winHeight - 25);
        glVertex2i(clamp_scaled_max + 4, winHeight - 25);
      }
		glEnd();
	}

  glEnable(GL_BLEND);
  glEnable(GL_LIGHTING);
  glEnable(GL_DEPTH_TEST);  
}

int renderer_get_grid_cell_size(void)
{
  return grid_cell_size;
}

void renderer_set_grid_cell_size(int cell_size)
{
  grid_cell_size = cell_size;
}

static void render_grid(void)
{
  int i, DIM;

  DIM = renderer_get_grid_cell_size();

  glColor4f(0.0f, 0.0f, 0.5f, 0.8f);
  glLineWidth(1.0f);

  glBegin(GL_LINES);
  for (i = 0; i < winWidth; i += DIM)
  {
    if (i)
    {
      glVertex3i(i, 0, 0);
      glVertex3i(i, winWidth, 0);
    }
  }

  for (i = 0; i < winHeight; i += DIM)
  {
    if (i)
    {
      glVertex3i(0, i, 0);
      glVertex3i(winHeight, i, 0);
    }
  }

  glEnd();
}

static void render_normal(void)
{
  int i, j, idx, DIM;
  double px, py;
	fftw_real  wn, hn;
  float normal_scale = 16.0f;

  DIM = fluids_get_dim();
  wn = (fftw_real)winWidth / (fftw_real)(DIM + 1);   // Grid cell width
	hn = (fftw_real)winHeight / (fftw_real)(DIM + 1);  // Grid cell height

  glDisable(GL_LIGHTING);
  glBegin(GL_LINES);
  glColor3f(1.0f, 0.2f, 0.0f);

  for (j = 0; j < DIM - 1; j++)
  {
    for (i = 0; i < DIM; i++)
    {
      idx = (j * DIM) + i;
      px = wn + (fftw_real)i * wn;
	    py = hn + (fftw_real)(j + 1) * hn;

      glVertex3f(px, py, height_array[idx]);
      glVertex3f(px + normal_array[idx].x * normal_scale,
                 py + normal_array[idx].y * normal_scale,
                 height_array[idx] + normal_array[idx].z * normal_scale);
    }
  }
  glEnable(GL_LIGHTING);
  glEnd();
}


static void render_normals(void)
{
  int i, j, idx, DIM;
	fftw_real  wn, hn;
  fftw_real *height;
  struct point *normal;

  DIM = fluids_get_dim();
  wn = (fftw_real)winWidth  / (fftw_real)(DIM + 1);  // Grid cell width
	hn = (fftw_real)winHeight / (fftw_real)(DIM + 1);  // Grid cell height

  height = heightplots_get_frame();
  normal = normals_get_frame();

  glDisable(GL_LIGHTING);
  glBegin(GL_LINES);
  glColor3f(0.2f, 1.0f, 0.0f);

  for (j = 0; j < DIM - 1; j++)
  {
    for (i = 0; i < DIM; i++)
    {
      int height_multiplier;
      double px, py, pz;

      height_multiplier = heightplots_get_height();
      idx = (j * DIM) + i;
      px = wn + (fftw_real)i * wn;
		  py = hn + (fftw_real)(j + 1) * hn;
      pz = height[idx] *height_multiplier;

      glVertex3f(px, py, pz);
      glVertex3f(px +(normal[idx].x *NORMAL_SCALE),
                 py +(normal[idx].y *NORMAL_SCALE),
                 pz +(normal[idx].z *NORMAL_SCALE));
    }
  }
  glEnable(GL_LIGHTING);
  glEnd();
}


void render_smoke(void)
{
  int        i, j, idx, DIM;
  double px, py, pz;
	fftw_real  wn, hn;
  struct color4f color;
  fftw_real *frame, *height;
  struct point *normal;

  DIM = fluids_get_dim();

  if (colormap_get_alpha() < 1.0f)
  {
    glEnable(GL_BLEND);
  }
  else
  {
    glDisable(GL_BLEND);
  }

  wn = (fftw_real)winWidth / (fftw_real)(DIM + 1);   // Grid cell width
  hn = (fftw_real)winHeight / (fftw_real)(DIM + 1);  // Grid cell height

  frame = smoke_get_frame();
  height = heightplots_get_frame();
  normal = normals_get_frame();

  glPolygonMode(GL_FRONT_AND_BACK, GL_FILL);

  for (j = 0; j < DIM - 1; j++)			//draw smoke
  {
	  glBegin(GL_TRIANGLE_STRIP);
	  i = 0;
	  idx = (j * DIM) + i;
	  px = wn + (fftw_real)i * wn;
	  py = hn + (fftw_real)j * hn;
    pz = height[idx] *heightplots_get_height();
		
    color = colormap_get_color(frame[idx]);
	  glColor4f(color.r, color.b, color.g, color.a);
//    glNormal3f(normal[idx].x, normal[idx].y, normal[idx].z);
	  glVertex3f(px, py, pz);

	  for (i = 0; i < DIM - 1; i++) 
	  {
		  idx = ((j + 1) * DIM) + i;
		  px = wn + (fftw_real)i * wn;
		  py = hn + (fftw_real)(j + 1) * hn;
      pz = height[idx] *heightplots_get_height();
      color = colormap_get_color(frame[idx]);
      glColor4f(color.r, color.g, color.b, color.a);      
  //    glNormal3f(normal[idx].x, normal[idx].y, normal[idx].z);
      glVertex3f(px, py, pz);
		  idx = (j * DIM) + (i + 1);
		  px = wn + (fftw_real)(i + 1) * wn;
		  py = hn + (fftw_real)j * hn;
      pz = height[idx] *heightplots_get_height();
		  color = colormap_get_color(frame[idx]);
      glColor4f(color.r, color.g, color.b, color.a);
   //   glNormal3f(normal[idx].x, normal[idx].y, normal[idx].z);
		  glVertex3f(px, py, pz);
	  }

	  idx = ((j + 1) * DIM) + (DIM - 1);
	  px = wn + (fftw_real)(DIM - 1) * wn;
	  py = hn + (fftw_real)(j + 1) * hn;
    pz = height[idx] *heightplots_get_height();
	  color = colormap_get_color(frame[idx]);
    glColor4f(color.r, color.g, color.b, color.a);
  //  glNormal3f(normal[idx].x, normal[idx].y, normal[idx].z);
	  glVertex3f(px, py, pz);
	  glEnd();
  }

  glDisable(GL_BLEND);
}


void render_grid_smoke(void)
{
  int i, j, idx, DIM;
  double px, py, pz;
	fftw_real  wn, hn;
  struct color4f color;
  fftw_real *frame, *height;

  DIM = fluids_get_dim();

  glDisable(GL_BLEND);

  wn = (fftw_real)winWidth / (fftw_real)(DIM + 1);   // Grid cell width
  hn = (fftw_real)winHeight / (fftw_real)(DIM + 1);  // Grid cell height

  frame = smoke_get_frame();
  height = heightplots_get_frame();

  for (j = 0; j < DIM - 1; j++)
  {
    glBegin(GL_LINE_STRIP);

    // horizontal lines
	  for (i = 0; i < DIM - 1; i++) 
	  {
		  // vertex 2
      idx = ((j + 1) * DIM) + i;
		  
      px = wn + (fftw_real)i * wn;
		  py = hn + (fftw_real)(j + 1) * hn;
      pz = height[idx] * heightplots_get_height();

		  color = colormap_get_color(frame[idx]);
      glColor4f(color.r, color.g, color.b, color.a);

      glVertex3f(px, py, pz);

		  // vertex 1
      idx = (j * DIM) + i;
		  
      px = wn + (fftw_real)i * wn;
		  py = hn + (fftw_real)j * hn;
      pz = height[idx] * heightplots_get_height();

      color = colormap_get_color(frame[idx]);
      glColor4f(color.r, color.g, color.b, color.a);

      glVertex3f(px, py, pz);

		  // vertex 3
      idx = (j * DIM) + (i + 1);
		  
      px = wn + (fftw_real)(i + 1) * wn;
		  py = hn + (fftw_real)j * hn;
      pz = height[idx] * heightplots_get_height();

		  color = colormap_get_color(frame[idx]);
      glColor4f(color.r, color.g, color.b, color.a);

      glVertex3f(px, py, pz);

	  }

    glEnd();
  }
}


static void render_glyph(GLUquadricObj *qobj, float x_value, float y_value, float i, float j)
{
  float x0, y0, z0, x1, y1, z1, x_dev, y_dev, size, length;
  float scale;
  double theta, in_prod;
  fftw_real wn, hn;
  int DIM;

  DIM = fluids_get_dim();
  scale = (float)((float)DIM / (float)fluids_get_var_dim()) / 6;
  size = quake_root((x_value * x_value * 20) + (y_value * y_value * 20)) * 5 * scale;

  wn = (fftw_real)winWidth / (fftw_real)(fluids_get_var_dim() + 1);  // Grid cell width
	hn = (fftw_real)winHeight / (fftw_real)(fluids_get_var_dim() + 1); // Grid cell heigh

  x0 = wn + (fftw_real)i * wn;
  y0 = hn + (fftw_real)j * hn;
  z0 = 16.0f;
  
  x1 = x0 + (x_value *1000)/4;
  y1 = y0 + (y_value *1000)/4;
  z1 = 0.0f;

  // inner product
  x_dev = x1 - x0;
  y_dev = y1 - y0;
  length = quake_root(x_dev * x_dev + y_dev * y_dev);
  
  x_dev = (length == 0) ? 0 : x_dev / length;
  y_dev = (length == 0) ? 1 : y_dev / length;
  
  in_prod = x_dev * 0 + y_dev * 1;
  theta = acos(in_prod) * (180/3.141592654);
  if (x1 > x0) { theta *= -1; }

  switch(glyphs_get_sort())
  {
    default:
    case GLYPH_LINES:
      glBegin(GL_LINES);
      glVertex3f(x0, y0, z0);
      glVertex3f(x1, y1, z0);
      glEnd();
      break;

    case GLYPH_TRIANGLES:
      if (size < 0.08f) return;

      glPushMatrix();
      
        glTranslatef(x0, y0, z0);
        glRotatef(theta, 0.0, 0.0, 1.0);
        glTranslatef(-x0, -y0, -z0);

        glBegin(GL_TRIANGLE_STRIP);
          glVertex3f(-10 * size + x0, -25 * size + y0, z0);
          glVertex3f( 0  * size + x0,  25 * size + y0, z0);
		      glVertex3f( 10 * size + x0, -25 * size + y0, z0);
        glEnd();

        glRotatef(-theta, 0.0, 0.0, 1.0);

      glPopMatrix();
      break;

    case GLYPH_CONES:
      x_dev = x1 - x0;
      y_dev = y1 - y0;
      length = quake_root(x_dev * x_dev + y_dev * y_dev) / 16;

      if (length < 0.5f) return;

      glPushMatrix();
      
        glTranslatef(x0, y0, z0);
        glRotatef(theta, 0.0, 0.0, 1.0);
        glRotatef(-90, 1.0, 0.0, 0.0);
        glScalef(length, length, length);
        gluCylinder(qobj, 10, 0.2, 8, 8, 8);

      glPopMatrix();
      break;

    case GLYPH_ARROWS:
    case GLYPH_QUAKE:
      if (size < 0.08f) return;

      glPushMatrix();
      
        glEnable(GL_TEXTURE_2D);
        
        glTranslatef(x0, y0, z0);
        glRotatef(theta, 0.0, 0.0, 1.0);
        glTranslatef(-x0, -y0, -z0);

        glBegin(GL_TRIANGLE_STRIP);
          glTexCoord2f(1, 1); glVertex3f( 20 * size + x0,  20 * size + y0, z0);
          glTexCoord2f(0, 1); glVertex3f(-20 * size + x0,  20 * size + y0, z0);
		      glTexCoord2f(1, 0); glVertex3f( 20 * size + x0, -20 * size + y0, z0);
          glTexCoord2f(0, 0); glVertex3f(-20 * size + x0, -20 * size + y0, z0);
        glEnd();

        glRotatef(-theta, 0.0, 0.0, 1.0);
        
        glDisable(GL_TEXTURE_2D);

      glPopMatrix();
      break;
  }
}

static void render_glyphs(void)
{
  int i, j, idx, DIM;
  float scale, idxcf, idxrf;
  struct color4f color;
  fftw_real *frame_color;
  struct fftw_real_xy frames_direction;
  GLUquadricObj *qobj;

  qobj = gluNewQuadric();
  gluQuadricDrawStyle(qobj, GLU_FILL);
  gluQuadricNormals(qobj, GLU_SMOOTH);

  DIM = fluids_get_dim();
  frame_color = glyphs_get_frame_color();
  frames_direction = glyphs_get_frames_direction();

  scale = (float)((float)DIM / (float)fluids_get_var_dim());
  for (i = 0; i < fluids_get_var_dim(); i++)
  {
    for (j = 0; j < fluids_get_var_dim(); j++)
    {
      idxcf = round(j * scale) * DIM;
      idxrf = round(i * scale);
      idx = idxcf + idxrf;

      color = glyphs_get_color(frame_color[idx]);
      glColor4f(color.r, color.b, color.g, color.a);

      render_glyph(qobj, frames_direction.x[idx], frames_direction.y[idx], i, j);
    }
  }

  gluDeleteQuadric(qobj);
}

struct point get_intersection(struct point pi, struct point pj, float fi, float fj, float t)
{
  struct point rp;
  rp.z = 8.0f;

  if (fi == min(fi, fj))
  {
    rp.x = (pi.x * (fj - t) + pj.x * (t - fi)) / (fj - fi);
    rp.y = (pi.y * (fj - t) + pj.y * (t - fi)) / (fj - fi);
  }
  else
  {
    rp.x = (pj.x * (fi - t) + pi.x * (t - fj)) / (fi - fj);
    rp.y = (pj.y * (fi - t) + pi.y * (t - fj)) / (fi - fj);
  }

  return rp;
}

static void render_isolines(void)
{
  int idx, count, i, j, state, DIM;
  float wn, hn, v0, v1, v2, v3, threshold, iso_scale;
  float x_offset, y_offset, z_pos = 8.0f;
  struct color4f color;
  struct point p0, p1, p2, p3, e1, e2;
  fftw_real *frame;

  DIM = fluids_get_dim();
  frame = isolines_get_frame();

  wn = (fftw_real)winWidth  / (fftw_real)(DIM + 1); // Grid cell width
  hn = (fftw_real)winHeight / (fftw_real)(DIM + 1); // Grid cell height

  if (isolines_get_nr())
  {
    iso_scale = fabs(isolines_get_threshold_min() - isolines_get_threshold_max()) / isolines_get_nr();
  }
  else
  {
    iso_scale = 0.0f;
  }

  for (count = 0; count < isolines_get_nr(); count++)
  {
    threshold = min(isolines_get_threshold_min(), isolines_get_threshold_max()) + count * iso_scale;

    glDisable(GL_LIGHTING);
    glLineWidth(2.0f);
    glBegin(GL_LINES);

    for (i = 0; i < DIM - 1; i++)
    {
      for (j = 0; j < DIM - 1; j++)
      {
        state = 0;
        idx = (j * DIM) + i;
        color = isolines_get_color(frame[idx]);
        glColor4f(color.r, color.g, color.b, color.a);

        v0 = frame[idx + DIM];
        v1 = frame[idx + 1 + DIM];
        v2 = frame[idx + 1];
        v3 = frame[idx];

        if (v0 >= threshold) { state += 1; }
        if (v1 >= threshold) { state += 2; }
        if (v2 >= threshold) { state += 4; }
        if (v3 >= threshold) { state += 8; }

        x_offset = wn + (fftw_real)i * wn;
        y_offset = (2 * hn) + (fftw_real)j * hn;

        p0.x = x_offset;
        p0.y = y_offset + hn;
        p0.z = z_pos;

        p1.x = x_offset + wn;
        p1.y = y_offset + hn;
        p1.z = z_pos;

        p2.x = x_offset + wn;
        p2.y = y_offset;
        p2.z = z_pos;

        p3.x = x_offset;
        p3.y = y_offset;
        p3.z = z_pos;

        switch(state)
        {
          default:
            state = 0;
          case 0:
          case 15:
            e1.x = e1.y = e1.z = 0.0f;
            e2.x = e2.y = e2.z = 0.0f;
            break;

          case 1:
          case 14:
            e1 = get_intersection(p0, p3, v0, v3, threshold);
            e2 = get_intersection(p0, p1, v0, v1, threshold);
            break;

          case 2:
          case 13:
            e1 = get_intersection(p1, p0, v1, v0, threshold);
            e2 = get_intersection(p1, p2, v1, v2, threshold);
            break;

          case 3:
          case 12:
            e1 = get_intersection(p0, p3, v0, v3, threshold);
            e2 = get_intersection(p1, p2, v1, v2, threshold);
            break;

          case 4:
          case 11:
            e1 = get_intersection(p2, p3, v2, v3, threshold);
            e2 = get_intersection(p2, p1, v2, v1, threshold);
            break;

          case 5:
          case 10:
            e1 = get_intersection(p0, p3, v0, v3, threshold);
            e2 = get_intersection(p0, p1, v0, v1, threshold);
            glVertex3i(e1.x, e1.y, e1.z);
            glVertex3i(e2.x, e2.y, e2.z);

            e1 = get_intersection(p2, p3, v2, v3, threshold);
            e2 = get_intersection(p2, p1, v2, v1, threshold);
            break;

          case 6:
          case 9:
            e1 = get_intersection(p2, p3, v2, v3, threshold);
            e2 = get_intersection(p1, p0, v1, v0, threshold);
            break;

          case 7:
          case 8:
            e1 = get_intersection(p2, p3, v2, v3, threshold);
            e2 = get_intersection(p0, p3, v0, v3, threshold);
            break;
        } // end switch state

        // draw line
        glVertex3i((GLint)e1.x, (GLint)e1.y, (GLint)e1.z);
        glVertex3i((GLint)e2.x, (GLint)e2.y, (GLint)e2.z);
      }
    } // end for count 
    glEnd();
  }
  glEnable(GL_LIGHTING);
}

static void render_streamlines(void)
{
  int i, j, k, l, idx, idx_x, idx_y, DIM;
  float v;
  struct point p;
  fftw_real *frame_history;
  fftw_real cell_x, cell_y;

  i = j = k = l = 0;
  DIM = fluids_get_dim();

  for (i = 0; i < get_cur_seedpoint(); i++)
  {
    p = get_seedpoint(i);
    
    glBegin(GL_LINE_STRIP);

    for (j = 0; j < HISTORY_SIZE; j++)
    {
      frame_history = streamlines_get_history(j);

      cell_x = (fftw_real)winWidth  / (fftw_real)(DIM + 1); // Grid cell width
	    cell_y = (fftw_real)winHeight / (fftw_real)(DIM + 1); // Grid cell heigh

      idx_x = round(p.x / cell_x);
      idx_y = round(p.y / cell_y);

      idx = (int)(cell_x * cell_y);

      v = (float)frame_history[idx];

      glVertex3f(p.x + v, winHeight - p.y, 0.0f - j * 10.0f);
    }
    glEnd();
  }
}

static int hisdex = 0;
 
void set_hisdex(int arg) {
  hisdex = arg;
}

static void render_flowvis(void)
{
	int        i, j, idx, DIM;
  double px,py;
	fftw_real  wn, hn;
  struct color4f color;
  fftw_real *history;

  
  DIM = fluids_get_dim();
  wn = (fftw_real)winWidth / (fftw_real)(DIM + 1);   // Grid cell width
	hn = (fftw_real)winHeight / (fftw_real)(DIM + 1);  // Grid cell height

  history = streamlines_get_history(hisdex);

  glPolygonMode(GL_FRONT_AND_BACK, GL_FILL);
  for (j = 0; j < DIM - 1; j++)			//draw smoke
  {
	  glBegin(GL_TRIANGLE_STRIP);

	  i = 0;
	  px = wn + (fftw_real)i * wn;
	  py = hn + (fftw_real)j * hn;
	  idx = (j * DIM) + i;
		
    color = flowvis_get_color((float)history[idx]);
    glColor4f(color.r, color.g, color.b, color.a);
	  glVertex2f((GLfloat)px, (GLfloat)py); // vertex 1
						
	  for (i = 0; i < DIM - 1; i++) 
	  {
		  px = wn + (fftw_real)i * wn;
		  py = hn + (fftw_real)(j + 1) * hn;
		  idx = ((j + 1) * DIM) + i;
      color = flowvis_get_color((float)history[idx]);
      glColor4f(color.r, color.g, color.b, color.a);
      glVertex2f((GLfloat)px, (GLfloat)py); // vertex 2

		  px = wn + (fftw_real)(i + 1) * wn;
		  py = hn + (fftw_real)j * hn;
		  idx = (j * DIM) + (i + 1);
		  color = flowvis_get_color((float)history[idx]);
      glColor4f(color.r, color.g, color.b, color.a);
		  glVertex2f((GLfloat)px, (GLfloat)py); // vertex 3
	  }

	  px = wn + (fftw_real)(DIM - 1) * wn;
	  py = hn + (fftw_real)(j + 1) * hn;
	  idx = ((j + 1) * DIM) + (DIM - 1);
	  color = flowvis_get_color((float)history[idx]);
    glColor4f(color.r, color.g, color.b, color.a);
	  glVertex2f((GLfloat)px, (GLfloat)py); // vertex 4

	  glEnd();
  }
}



void renderer_init_gl(void)
{
  float LightAmbient[]  = { 0.10f, 0.10f, 0.10f, 1.0f }; // Ambient light values
  float LightDiffuse[]  = { 0.80f, 0.80f, 0.80f, 1.0f }; // Diffuse light values
  float LightPosition[] = { x_pos, y_pos, -500, 1.0f }; // Position of the light source

  glClearColor(0.0f, 0.0f, 0.0f, 0.0f);
  glHint(GL_PERSPECTIVE_CORRECTION_HINT, GL_NICEST);
  glHint(GL_POINT_SMOOTH_HINT, GL_NICEST);

  // blending
  glEnable(GL_BLEND);
  glBlendFunc(GL_SRC_ALPHA, GL_ONE);

  // depth testing
  glEnable(GL_DEPTH_TEST);
  glClearDepth(1.0f);
  glDepthFunc(GL_LESS);
  
  // lighting
  glShadeModel(GL_SMOOTH);
  glEnable(GL_LIGHTING);
  glEnable(GL_LIGHT0);
  glEnable(GL_COLOR_MATERIAL);

  glMaterialf(GL_FRONT_AND_BACK, GL_SHININESS, 9.8f);

  glLightfv(GL_LIGHT0, GL_AMBIENT,  LightAmbient);
  glLightfv(GL_LIGHT0, GL_DIFFUSE,  LightDiffuse);
  glLightfv(GL_LIGHT0, GL_POSITION, LightPosition);
}


void renderer_set_render_grid(int render_grid)
{
  renderer_render_grid = render_grid;
}

int renderer_get_render_grid(void)
{
  return renderer_render_grid;
}

void renderer_set_zoomspeed(int zoomspeed)
{
  renderer_zoomspeed = zoomspeed;
}

int renderer_get_zoomspeed(void)
{
  return renderer_zoomspeed;
}

void renderer_zoom_in(void)
{
  if ((z_pos +renderer_zoomspeed) < MIN_ZOOM) {
    z_pos += renderer_zoomspeed;
  }
}

void renderer_zoom_out(void)
{
  if ((z_pos -renderer_zoomspeed) > MAX_ZOOM) {
    z_pos -= renderer_zoomspeed;
  }
}

void renderer_reset_zoom(void)
{
  z_pos =  DEFAULT_ZOOM;
}

void renderer_reset_pos(void)
{
  x_pos = DEFAULT_X_POS;
  y_pos = DEFAULT_Y_POS;
}

void renderer_reset_yaw(void)
{
  x_rot = 0.0f;
}

void renderer_reset_pitch(void)
{
  y_rot = 0.0f;
}

void renderer_reset_roll(void)
{
  z_rot = 0.0f;
}


//visualize: This is the main visualization function
void visualize(struct vis_data_arrays *vis_data)
{
  // Rotate field
  glPushMatrix();

    glLoadIdentity();
    glTranslatef(x_pos, y_pos, z_pos);
    glRotatef(x_rot, 1.0f, 0.0f, 0.0f);
    glRotatef(y_rot, 0.0f, 1.0f, 0.0f);
    glRotatef(z_rot, 0.0f, 0.0f, 1.0f);
  
    glTranslatef((GLfloat)(-winWidth / 2), (GLfloat)(-winHeight / 2), 0.0f);

    if (renderer_render_grid) {
      render_grid();
    }
    if (smoke_get_render())
    {
      if (smoke_get_sort() == SMOKE_FLUID)
      {
	      render_smoke();
      }
      else // SMOKE_GRID
      {
        render_grid_smoke();
      }
    }

    if (isolines_get_render()) {
      render_isolines();
    }

    if (glyphs_get_render()) {
      render_glyphs();
    }

    render_seedpoints();

    if (streamlines_get_render()) {
      render_streamlines();
    }

    if (flowvis_get_render()) {
      render_flowvis();
    }

    if (normals_get_render()) {
      render_normal();
    }

    if (normals_get_render2()) {
      render_normals();
    }

    render_streamlines();
    //create an orthographic projection matrix for rendering the legendbar
    glMatrixMode(GL_PROJECTION);
    glPushMatrix();

      glLoadIdentity();
      
      //create an orthographic viewing volume with the same dimensions as the 
      //window in front of the previous viewing volume
      //from znear = 0.0 to zfar = 1.0
      glOrtho(0.0, (GLfloat)winWidth, 0.0, (GLfloat)winHeight, 0.0, 1.0);
      glMatrixMode(GL_MODELVIEW);
      glPushMatrix(); //save the current modelview matrix

        glLoadIdentity();

        //now render the legendbar
        render_legend();

      //restore the previous modelview matrix
      glPopMatrix();

      //restore the previous projection matrix
      glMatrixMode(GL_PROJECTION);

    glPopMatrix();
    glMatrixMode(GL_MODELVIEW);

  glPopMatrix();
}