Visualization with OpenGL hands-on

Tutorial aims to introduce visualization techniques with modern Open Graphics Language (OpenGL) approaches standardized with version 3.x+. OpenGL Shading Language (GLSL) is used for that without tendency to introduce photo-realismas output but rather useful colors for scientific data exploration.

Running this tutorial on Linux desktop one requires at least the OpenGL 2.0 graphics with the GLSL 1.1 and supporting libraries GL, GLU, GLUT, GLEW. This can be verified with the following commands:

$ glxinfo |grep OpenGL.*version
OpenGL version string: 2.1 Mesa 8.0.5
OpenGL shading language version string: 1.20
$ ls /usr/include/GL/{glut.h,glew.h,gl.h,glu.h}
/usr/include/GL/glew.h  /usr/include/GL/glu.h
/usr/include/GL/gl.h    /usr/include/GL/glut.h

Legacy OpenGL

Create the following first.c using your favorite editor.

#include <GL/glut.h>

void display()
{
  glClear(GL_COLOR_BUFFER_BIT);
  glColor3f(1.0, 0.4, 1.0);
  glBegin(GL_LINES);
    glVertex2f(0.1, 0.1);
    glVertex3f(0.8, 0.8, 1.0);
  glEnd();
  glutSwapBuffers();
}

int main(int argc, char *argv[])
{
  glutInit(&argc,argv);
  glutInitDisplayMode(GLUT_DOUBLE);
  glutCreateWindow("first.c GL code");
  glutDisplayFunc(display);
  glutMainLoop();
  return 0;
}

Create Makefile to build your program.

CFLAGS=-Wall -g
LDFLAGS=-lGL -lGLU -lglut -lGLEW

ALL=first
default: $(ALL)

first : first.o

clean:
      rm -rf  *.o *~ [!m]*.obj core* $(ALL)

Beware that Makefile is TAB aware. So the last line should contain TAB indentation and not spacing.

Make and run the program with

make
./first

Try the same program in Python

from OpenGL.GLUT import *
from OpenGL.GL import *
import sys

def display():
    glClear(GL_COLOR_BUFFER_BIT)
    glColor3f(1.0, 0.4, 1.0)
    glBegin(GL_LINES)
    glVertex2f(0.1, 0.1)
    glVertex3f(0.8, 0.8, 1.0)
    glEnd()
    glutSwapBuffers()

if __name__ == "__main__":
    glutInit(sys.argv)
    glutInitDisplayMode(GLUT_DOUBLE)
    glutCreateWindow("first.py GL code")
    glutDisplayFunc(display)
    glutMainLoop()

and run it with

python first.py

Exercises #1:

1. Add RGB color to vertices with glColor3f(0.0, 0.4, 1.0);.
2. Replace single line drawing in display() with the following snippet

      GLfloat vertices[][2] = {
       { -0.90, -0.90 }, // Triangle 1
       {  0.85, -0.90 },
       { -0.90,  0.85 },
       {  0.90, -0.85 }, // Triangle 2
       {  0.90,  0.90 },
       { -0.85,  0.90 }
      };
   

   and try to draw two wireframe triangles in a loop. Change primitive to GL_LINE_LOOP.
3. Draw two primitives with GL_TRIANGLES.
4. Add different color to each vertex.

    GLfloat color[][3] = {
     {1, 0, 0}, {0, 1, 0}, {0, 0, 1},
     {1, 1, 0}, {0, 1, 1}, {1, 0, 1}};
   

5. Replace loop with the following

     glVertexPointer(2, GL_FLOAT, 0, &vertices[0][0]);
     glEnableClientState(GL_VERTEX_ARRAY);
     glDrawArrays(GL_TRIANGLES, 0, 6);
     glDisableClientState(GL_VERTEX_ARRAY);
   

   How can we add color to vertices? See ​glColorPointer and ​glEnableClientState.
6. Change background to glClearColor(0.9,1,1,1.0); and suggest initial window in main()

   glutInitWindowSize(512, 512);
   glutInitWindowPosition((glutGet(GLUT_SCREEN_WIDTH)-512)/2,
                          (glutGet(GLUT_SCREEN_HEIGHT)-512)/2);
   

7. Add ​keyboard event to quit the program when pressing ESCape key with keycode 27 by adding callback function

    void keyboard(unsigned char key, int x, int y) 
    {
      if (key == 27) 
             exit(0);
    }
   

   and registering event within main() by glutKeyboardFunc(keyboard);. Some prefer key == 'q', though.

Modern OpenGL

We extend previous exercise with example that introduces OpenGL 3.x techniques:

• OpenGL Shading Language (GLSL 1.2) where simple vertex and fragment shader are required.
• Vertex Aray Objects (VAOs) and vertex buffers (VBOs) stored in GPU.

Create triangle.c and update Makefile with new target

#include <stdio.h>
#include <stdlib.h>
#include <GL/glew.h>
#include <GL/glut.h>

GLuint program;
GLuint vbo_vertices;
GLint attribute_coord2d;

static const GLchar * vertex_shader[] = {
  "attribute vec2 coord2d;" // input vertex position
  "void main()"
  "{"
  "  gl_Position = gl_ModelViewProjectionMatrix*vec4(coord2d, 0.0, 1.0);"
  "}"
};
static const GLchar * fragment_shader[] =
  {"void main()"
   "{"
   " gl_FragColor = vec4(0.4,0.4,0.8,1.0);"
   "}"
  };

void create_shaders()
{
  GLuint v, f;

  v = glCreateShader(GL_VERTEX_SHADER);
  f = glCreateShader(GL_FRAGMENT_SHADER);
  glShaderSource(v, 1, vertex_shader, NULL);
  glShaderSource(f, 1, fragment_shader, NULL);
  glCompileShader(v);
  glCompileShader(f);
  program = glCreateProgram();
  glAttachShader(program, f);
  glAttachShader(program, v);
  glLinkProgram(program);
  glUseProgram(program);

  attribute_coord2d = glGetAttribLocation(program, "coord2d");
  if (attribute_coord2d == -1) {
    fprintf(stderr, "Could not bind attribute coord2d\n");
  }
  glEnableVertexAttribArray(attribute_coord2d);
}

void send_buffers_to_GPU(void)
{
  GLuint vertex_array_object;
  glGenVertexArrays(1, &vertex_array_object);
  glBindVertexArray(vertex_array_object);
  
  GLfloat vertices[][2] = {
    { -0.90, -0.90 }, // Triangle 1
    {  0.85, -0.90 },
    { -0.90,  0.85 },
    {  0.90, -0.85 }, // Triangle 2
    {  0.90,  0.90 },
    { -0.85,  0.90 }
  };

  glGenBuffers(1, &vbo_vertices);
  glBindBuffer( GL_ARRAY_BUFFER, vbo_vertices);
  glBufferData( GL_ARRAY_BUFFER, sizeof(vertices), vertices, GL_STATIC_DRAW);
}


void display(void)
{
  glClear(GL_COLOR_BUFFER_BIT);
  glBindBuffer(GL_ARRAY_BUFFER, vbo_vertices);
  glVertexAttribPointer(attribute_coord2d, 2, GL_FLOAT, GL_FALSE, 0, NULL);
  glDrawArrays(GL_TRIANGLES, 0, 6); // Draw 6 vertices
  glutSwapBuffers();
}

int main(int argc, char **argv)
{
  glutInit(&argc, argv);
  glutInitDisplayMode(GLUT_RGBA | GLUT_DOUBLE);
  glutCreateWindow("GLSL Intro");
  glutDisplayFunc(display);
  glewInit();
  if (!glewIsSupported("GL_VERSION_2_0"))
   {
     printf("GLSL not supported\n");
     exit(EXIT_FAILURE);
   }
  glClearColor(0.9,1.0,1.0,1.0);
  send_buffers_to_GPU();
  create_shaders();
  glutMainLoop();
  return EXIT_SUCCESS;
}

Exercises #2

1. To be able to continue and not get lost introduce shader compiler logs in case of compilation errors by adding the following code into create_shaders() right at after vertex shader compilation:

   GLint compiled;
   glGetShaderiv(v, GL_COMPILE_STATUS, &compiled );
   if ( !compiled ) {
     GLsizei  maxLength, length;
     glGetShaderiv( v, GL_INFO_LOG_LENGTH, &maxLength );
     GLchar* log = malloc(sizeof(GLchar)*(maxLength+1));
     glGetShaderInfoLog(v,  maxLength, &length, log);
     printf("Vertex Shader compilation failed: %s\n", log);
     free(log);
   }
   

   Do not forget to repeat the same thing for fragment shader. Add linker debugging

   GLint linked;
   glGetProgramiv(program, GL_LINK_STATUS, &linked );
   if ( !linked ) {
     GLsizei len;
     glGetProgramiv(program, GL_INFO_LOG_LENGTH, &len );
     GLchar* log = malloc(sizeof(GLchar)*(len+1));
     glGetProgramInfoLog(program, len, &len, log );
     printf("Shader linking failed: %s\n", log);
     free(log);
   }
   

   Create some error to verify if it works. For general (core) OpenGL errors we can use the following glcheck() utility at suspicious places.

   #define glcheck() {GLenum s; if ((s=glGetError()) != GL_NO_ERROR) \
                   fprintf (stderr, "OpenGL Error: %s at line %d\n", \
                            gluErrorString(s),  __LINE__);}
   

2. Copy triangle.c into temperature.c and introduce vertex temperature with additional array and buffer at the end of send_buffers_to_GPU()

   GLfloat vertex_temperature[] = {0, 1, 0.2, 0.1, 0.5, 0.9};
   glGenBuffers(1, &vbo_temperature);
   glBindBuffer(GL_ARRAY_BUFFER, vbo_temperature);
   glBufferData(GL_ARRAY_BUFFER, sizeof(vertex_temperature),
              vertex_temperature, GL_STATIC_DRAW);
   

   and adding corresponding global attribute and VBOs IDs at the top of temperature.c so that global section reads:

   GLuint program;
   GLuint vbo_vertices;
   GLuint vbo_temperature;
   GLint attribute_coord2d;
   GLint attribute_temperature;
   

   Replace shaders with

   static const GLchar * vertex_shader[] = {
     ""
     "attribute float temperature;" // custom variable along with vertex position
     "varying float  t;" // communicate between the vertex and the fragment shader
     "void main()"
     "{"
     "  t = temperature;"
     "  gl_Position = gl_ModelViewProjectionMatrix * gl_Vertex;"
     "}"
   };
   static const GLchar * fragment_shader[] = {
     "vec3 Cool = vec3(0, 0, 1);" // Red
     "vec3 Hot  = vec3(1, 0, 0);" // Blue
     "varying float t;" // Interpolated by fragment
     "void main()"
     "{"
     "  vec3 color = mix(Cool, Hot, t);"  // use the built-in mix() function
     "  gl_FragColor = vec4(color, 1.0);" // append alpha channel
     "}"
   };
   

   Bind temperature buffer and specify data layout within display() just before glDrawElements() with

   glBindBuffer(GL_ARRAY_BUFFER, vbo_temperature);
   glVertexAttribPointer(attribute_temperature, 1, GL_FLOAT, GL_FALSE, 0, NULL);
   

   What happens if we don't enable temperature vertex array? Confer temperature.c​ attached if having troubles.
3. Add additional custom vertex array for the pressure. Change the temperature array to have values in Celsius for water boiling range [0-100]°C. Pressure should be in the range of [0-1] MPa. Scaling to color range [0-1] should be done in shaders. Toggle between them with the keyboard event by using the keys 'p' and 't' that glEnableVertexAttribArray() and glDisableVertexAttribArray() corresponding vertex attribute arrays.

Interactivity

Assemble the following Utah teapot model and attached virtual trackball.h​ and trackball.c​ sources from SGI.

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

#include <GL/glew.h>
#include <GL/glut.h>
#include "trackball.h"

GLuint program;

static const GLchar * vertex_shader[] ={"\
varying vec3 normal, lightDir;\
uniform mat4 RotationMatrix;\
void main()\
{          \
  lightDir=normalize(vec3(gl_LightSource[0].position));\
  normal=normalize(gl_NormalMatrix*gl_Normal);\
  gl_Position = gl_ProjectionMatrix * \
  RotationMatrix*gl_ModelViewMatrix*gl_Vertex;\
}"};

static const GLchar * fragment_shader[] ={"\
/* simple toon fragment shader */\
/* www.lighthouse3d.com        */\
\
varying vec3 normal, lightDir;\
\
void main()\
{\
        float intensity;\
        vec3 n;\
        vec4 color;\
\
        n = normalize(normal);\
        intensity = max(dot(lightDir,n),0.0);\
        if (intensity > 0.98)\
                color = vec4(0.8,0.8,0.8,1.0);\
        else if (intensity > 0.5)\
                color = vec4(0.4,0.4,0.8,1.0);\
        else if (intensity > 0.25)\
                color = vec4(0.2,0.2,0.4,1.0);\
        else\
                color = vec4(0.1,0.1,0.1,1.0);\
        gl_FragColor = color;\
}"};

void create_shaders()
{
  GLuint v, f;

  v = glCreateShader(GL_VERTEX_SHADER);
  f = glCreateShader(GL_FRAGMENT_SHADER);
  glShaderSource(v, 1, vertex_shader, NULL);
  glShaderSource(f, 1, fragment_shader, NULL);
  glCompileShader(v);
  GLint compiled;
  glGetShaderiv(v, GL_COMPILE_STATUS, &compiled );
  if ( !compiled ) {
    GLsizei  maxLength, length;
    glGetShaderiv( v, GL_INFO_LOG_LENGTH, &maxLength );
    GLchar* log = malloc(sizeof(GLchar)*(maxLength+1));
    glGetShaderInfoLog(v,  maxLength, &length, log);
    printf("Vertex Shader compilation failed: %s\n", log);
    free(log);
  }
  glCompileShader(f);
  glGetShaderiv(f, GL_COMPILE_STATUS, &compiled );
  if ( !compiled ) {
    GLsizei  maxLength, length;
    glGetShaderiv( f, GL_INFO_LOG_LENGTH, &maxLength );
    GLchar* log = malloc(sizeof(GLchar)*(maxLength+1));
    glGetShaderInfoLog(f,  maxLength, &length, log);
    printf("Fragment Shader compilation failed: %s\n", log);
    free(log);
  }
  program = glCreateProgram();
  glAttachShader(program, f);
  glAttachShader(program, v);
  glLinkProgram(program);
  GLint linked;
  glGetProgramiv(program, GL_LINK_STATUS, &linked );
  if ( !linked ) {
    GLsizei len;
    glGetProgramiv(program, GL_INFO_LOG_LENGTH, &len );
    GLchar* log = malloc(sizeof(GLchar)*(len+1));
    glGetProgramInfoLog(program, len, &len, log );
    printf("Shader linking failed: %s\n", log);
    free(log);
  }
  glUseProgram(program);
}


float lpos[4] = {1, 0.5, 1, 0};
GLfloat m[4][4]; // modelview rotation matrix
float last[4], cur[4]; // rotation tracking quaternions 
int width, height, beginx, beginy;
float p1x, p1y, p2x, p2y;

void display(void) {
  GLuint location = glGetUniformLocation(program, "RotationMatrix");
  build_rotmatrix(m, cur);  
  glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
  glLightfv(GL_LIGHT0, GL_POSITION, lpos);
  if( location >= 0 )
    glUniformMatrix4fv(location, 1, GL_FALSE, &m[0][0]);
  glutSolidTeapot(0.6);
  glutSwapBuffers();
}

void reshape (int w, int h)
{
  double l = 1;
  width=w;  height=h;
  glViewport (0, 0, w, h);
  glMatrixMode (GL_PROJECTION);
  glLoadIdentity();
  glOrtho(-l, l, -l, l, -l, l);
  glMatrixMode(GL_MODELVIEW);
  glLoadIdentity();
}

void keys(unsigned char key, int x, int y)
{
   if (key == 27 || key == 'q') 
         exit(0);
}

void mouse(int button,int state, int x, int y)   
{
  beginx = x;
  beginy = y;
}

void motion(int x,int y)   
{
  p1x = (2.0*beginx - width)/width;
  p1y = (height - 2.0*beginy)/height;
  p2x = (2.0 * x - width) / width;
  p2y = (height - 2.0 * y) / height;
  trackball(last, p1x, p1y, p2x, p2y);   
  add_quats(last, cur, cur);   
  beginx = x;
  beginy = y;
  glutPostRedisplay();   
}

int main(int argc, char **argv)
{
  glutInit(&argc, argv);
  glutInitDisplayMode(GLUT_DEPTH | GLUT_DOUBLE | GLUT_RGBA);
  glutInitWindowSize(512, 512);
  glutInitWindowPosition((glutGet(GLUT_SCREEN_WIDTH)-512)/2,
                         (glutGet(GLUT_SCREEN_HEIGHT)-512)/2);
  glutCreateWindow("Use mouse to rotate");
  
  trackball(cur, 0.0, 0.0, 0.0, 0.0);

  glutDisplayFunc(display);
  glutReshapeFunc(reshape);
  glutMouseFunc(mouse);
  glutMotionFunc(motion);
  glutKeyboardFunc(keys);

  glEnable(GL_DEPTH_TEST);
  glClearColor(1.0,1.0,1.0,1.0);
  glewInit();
  if (!glewIsSupported("GL_VERSION_2_0"))
   {
     printf("GLSL not supported\n");
     exit(EXIT_FAILURE);
   }
  create_shaders();
  glutMainLoop();
  return EXIT_SUCCESS;
}

To build two sources we add the following line to Makefile:

teapot: teapot.o trackball.o

Exercises #3

1. OK, it rotates. But how come the light rotates with the teapot? I'm pretty sure that the light is not rotated while the teapot vertices are. Answer: Take a look into the normal. You need to rotate the vertex normal in the vertex_shader[] too! With a code like:

   vec4 n = RotationMatrix*vec4(gl_NormalMatrix*gl_Normal, 1);\
   normal = normalize(n.xyz); \
   

   
2. Introduce zoom in/out functionality of the viewer by adding mouse wheel events to the end of mouse()

   if (button == 3 && state == GLUT_DOWN)
     { zoom *= 1.1; glutPostRedisplay(); }
   else if (button == 4 && state == GLUT_DOWN)
     { zoom /= 1.1; glutPostRedisplay(); }
   

   and introduction of global variable float zoom = 1.0; that is communicated to GPU by additional uniform float Zoom; in the vertex_shader[]. Last line is replaced by

     gl_Position = gl_ProjectionMatrix * RotationMatrix \
      * gl_ModelViewMatrix*vec4(Zoom*gl_Vertex.xyz, 1.0); \
   

   In the display() we send zoom to GPU before drawing the glutSolidTeapot() by adding

    location = glGetUniformLocation(p, "Zoom");
    if (location >= 0)  glUniform1f(location, zoom);
   

3. Simplify the cartoon shader to ​Gouraud shader by using just greyscale gl_FragColor = vec4(intensity); or copper mix

    vec4 copper_ambient = vec4(0.191250, 0.073500, 0.022500, 1.000000); \
    vec4 copper_diffuse = vec4(0.703800, 0.270480, 0.082800, 1.000000); \
    gl_FragColor = copper_ambient + intensity*copper_diffuse;\
   

Reading Objects

Sometimes we hit limitations of the visualisation tools for the data that we want to visualize. For example motorBike.obj​ from ​OpenFOAM contains object groups that we want to show colored separately and not as whole. Neither ​VisIt and ​ParaView can read ​Wavefront OBJ file with group separation. We are forced to convert motorBike.obj​ into bunch of files and read them one by one. The following wavefront.c converts motorBike.obj​ into 67 files. Try to open them in ​VisIt and ​ParaView. Note that we need to compensate vertex counting that starts with 1 and not with 0.

#include <stdio.h>

#define MaxVertices 400000
#define MaxFaces    400000 
#define MaxGroups   100

float       vertex[MaxVertices*3];
unsigned int  face[MaxFaces*3];
char    group_name[MaxGroups][80];
int     start_face[MaxGroups];

int vertices = 0;
int faces    = 0;
int groups   = 0;

void read_wavefront(const char *filename)
{
  char line[80];
  FILE *f = fopen(filename, "r");
  while(fgets(line, sizeof(line), f))
    switch(line[0])
      {
      case 'v':
        sscanf(&line[1],  "%f %f %f", &vertex[vertices+3],
               &vertex[vertices*3+1], &vertex[vertices*3+2]);
        ++vertices;
        break;
      case 'g':
        sscanf(&line[1], "%s", group_name[groups]);
        start_face[groups++] = faces;
        break;
      case 'f':
        sscanf(&line[1],  "%d %d %d", &face[faces*3],
               &face[faces*3+1], &face[faces*3+2]);
        --face[faces*3]; --face[faces*3+1];
        --face[faces*3+2]; ++faces; 
        break;
      }
  fclose(f);
  start_face[groups] = faces;
  printf("Read %d vertices and %d faces within %d groups from %s\n",
         vertices, faces, groups, filename);
}

void write_wavefront(int group_number)
{
  int i = 0; char n[80], *p = group_name[group_number];
  while (*p != '%' && *p != '\0') n[i++] = *p++; // remove % from name
  n[i++] = '.'; n[i++] = 'o'; n[i++] = 'b'; n[i++] = 'j'; n[i] = '\0';
  FILE *f = fopen(n, "w"); fprintf(f, "# Wavefront OBJ file\n");
  for (i = 0; i < vertices; i++)
    fprintf(f, "v %g %g %g\n", vertex[i*3], vertex[i*3+1], vertex[i*3+2]);
  fprintf(f, "g %s\n", group_name[group_number]);
  for (i = start_face[group_number]; i < start_face[group_number+1]; ++i)
    fprintf(f, "f %d %d %d\n", face[i*3]+1, face[i*3+1]+1, face[i*3+2]+1);
  fclose(f);
}

int main(int argc, char **argv)
{
  int i;
  read_wavefront("motorBike.obj");
  for(i = 0; i < groups; i++) write_wavefront(i);
  return 0;
}

Exercises #4

1. Insert teapot.c interactivity example into wavefront.c above and save it as motorbike.c. Verify that there are no compile problems and that the main() contains read_wavefront("motorBike.obj");. Disable lengthy wavefront saving in main().
2. Instead of drawing of the teapot in display() we will draw single vertex array object as a point cloud with

    //glutSolidTeapot(0.6);
    glDrawArrays(GL_POINTS, 0, vertices);
   

   To be able to draw this, we need to prepare Vertex Array Object and buffer with init(); that is called from main() just before glutMainLoop();

   void init()
   {
     GLuint VAO[1];
     GLuint Buffer[1];
     glGenVertexArrays(1, VAO);
     glBindVertexArray(VAO[0]);
     glGenBuffers(1, Buffer);
     glBindBuffer(GL_ARRAY_BUFFER, Buffer[0]);
     glBufferData(GL_ARRAY_BUFFER, vertices*3*sizeof(GLfloat),
                  vertex, GL_STATIC_DRAW);
     glVertexAttribPointer(0, 3, GL_FLOAT, GL_FALSE, 0, NULL);
     glEnableVertexAttribArray(0);
   }
   

3. We see that rotation of the motorbike around the front wheel is not really nice. To compensate we translate all points in the vertex_shader[] by adding vec3(-0.75, 0, -0.7); to every vertex in world coordinates and thus moving motor bike to origin. Last part of the vertex_shader[]now reads:

   vec3 position = gl_Vertex.xyz + vec3(-0.75, 0, -0.7); \
   gl_Position = gl_ProjectionMatrix * RotationMatrix \
    * gl_ModelViewMatrix*vec4(Zoom*position, 1.0); \