Line Integral Convolution

Brian Cabral and Casey Leedom have developed a technique for vector field visualization known as line integral convolution.[8] The technique takes an input vector field and an input image. For each location p in the input vector field, a parametric curve P(p, s) is generated which passes through the location and follows the vector field for some distance in either direction. To create an output pixel F'(p), a weighted sum of the values of the input image F along the curve is computed. The weighting function is k(x). Thus the continuous form of the equation is:

To discretize the equation, we use values along the curve P(p,s):

The computation of the output values of this equation may accelerated using OpenGL. We use a mesh texture mapped with the input image to create the output image. The mesh is redrawn l times. At each step, we advect the texture coordinates and accumulate the results. Advection applies a mapping defined by the vector field to the input points. A simple advection implementation moves each point by a fixed amount in the direction of the vector flow at the point. Advection has been well-studied, and many more complicated algorithms exist.

An implementation of the algorithm uses the following variables:

• int l: Number of steps
• GLfloat h[0..(l-1)]: Kernel weights
• GLfloat hNormalize: Normalization factor ()
• GLfloat gridW, gridH: Size of the grid
• GLfloat *grid[2]: Grid of texture coordinates.

and the functions:

• advect_grid(GLfloat s): Advect grid by s, which may be positive or negative.

We begin by initializing the grid points:

void init(void)
{
  int x, y;
  
  for (y = 0; y < gridH; y++) {
    for (x = 0; x < gridW; x++) {
      grid[y*gridW + x][0] = x;
      grid[y*gridW + x][1] = y;
    }
  }
}

The texture image is then downloaded and bound. In the draw routine we call:

void lic(void)
{
  int x, y;
  int i;
  
  advect_grid(-l/2);

  glClear(GL_COLOR_BUFFER_BIT | GL_ACCUM_BUFFER_BIT);

  /* scale texture coordinates */
  glPushAttrib(GL_TRANSFORM_BIT);
  glMatrixMode(GL_TEXTURE);
  glPushMatrix();
  glScalef(1.0/(gridW-1), 1.0/(gridH-1), 1);

  for (i = 0; i < l; i++) {
    glEnable(GL_TEXTURE_2D);
    for (y = 0; y < gridH-1; y++) {
      glBegin(GL_QUAD_STRIP);
      for (x = 0; x < gridW-1; x++) {
        glTexCoord2fv(grid[y*gridW + x]);
        glVertex2i(x, y);
        glTexCoord2fv(grid[y*gridW + x+1]);
        glVertex2i(x+1, y);
        glTexCoord2fv(grid[(y+1)*gridW + x]);
        glVertex2i(x, y+1);
        glTexCoord2fv(grid[(y+1)*gridW + x+1]);
        glVertex2i(x+1, y+1);
      }
      glEnd();
    }
    glDisable(GL_TEXTURE_2D);
    glAccum(GL_ACCUM, h[i]);

    advect_grid(1);
  }

  glAccum(GL_RETURN, hNormalize);

  glPopMatrix();
  glPopAttrib();
}

In the lic routine, we first clear the color and accumulation buffers. Next, we modify the texture matrix such that a texture coordinate of (gridW, gridH) will map to the upper right corner of the input texture.

Upon each iteration of the loop, we draw the grid using the array of texture coordinates (vertex arrays could provide a more efficient implementation). Then, we accumulate the results, weighting by the kernel array entry. Next, we call advect_grid to update the texture coordinate array. At the end of the routine, we return the results and normalize by the sum of the kernel weights.

Upon implementation, several difficulties may present themselves. First, implementing advect_grid well is non-trivial (but well-studied). Second, here we have used a static grid to draw the field. This approach will probably lead to artifacts when drawing high-frequency fields or unnecessary inefficiency when drawing low-frequency fields. A better approach would subdivide the grid based on the behavior of the vector field. Also, the user may find that the results of the accumulation operation go outside the range [-1..1] if care is not taken when choosing the kernel and normalization values. Finally, dealing with the three different coordinate spaces (vector field, grid, and texture image) can become complicated.