git.hungrycats.org Git - xscreensaver/blob - hacks/glx/marching.c

   1 /* xscreensaver, Copyright (c) 2002 Jamie Zawinski <jwz@jwz.org>
   2  * Utility functions to create "marching cubes" meshes from 3d fields.
   3  *
   4  * Permission to use, copy, modify, distribute, and sell this software and its
   5  * documentation for any purpose is hereby granted without fee, provided that
   6  * the above copyright notice appear in all copies and that both that
   7  * copyright notice and this permission notice appear in supporting
   8  * documentation.  No representations are made about the suitability of this
   9  * software for any purpose.  It is provided "as is" without express or
  10  * implied warranty.
  11  *
  12  * Marching cubes implementation by Paul Bourke <pbourke@swin.edu.au>
  13  * http://astronomy.swin.edu.au/~pbourke/modelling/polygonise/
  14  */
  15
  16 #ifdef HAVE_CONFIG_H
  17 # include "config.h"
  18 #endif
  19
  20 #include <stdlib.h>
  21 #include <stdio.h>
  22 #include <math.h>
  23
  24 #ifdef HAVE_COCOA
  25 # include <OpenGL/gl.h>
  26 #else
  27 # include <GL/gl.h>
  28 #endif
  29
  30 #include "marching.h"
  31 #include "normals.h"
  32
  33 extern char *progname;
  34
  35 #undef ABS
  36 #define ABS(x) ((x)<0?(-(x)):(x))
  37
  38 typedef struct {
  39    XYZ p[3];
  40 } TRIANGLE;
  41
  42 typedef struct {
  43    XYZ p[8];
  44    double val[8];
  45 } GRIDCELL;
  46
  47
  48 /* Indexing convention:
  49
  50              Vertices:                    Edges:
  51
  52           4  ______________ 5           ______________
  53            /|            /|           /|     4      /|
  54           / |         6 / |       7  / |8        5 / |
  55       7  /_____________/  |        /______________/  | 9
  56         |   |         |   |        |   |   6     |   |
  57         | 0 |_________|___| 1      |   |_________|10_|
  58         |  /          |  /      11 | 3/     0    |  /
  59         | /           | /          | /           | / 1
  60       3 |/____________|/ 2         |/____________|/
  61                                           2
  62  */
  63
  64 static const int edgeTable[256] = {
  65   0x0  , 0x109, 0x203, 0x30a, 0x406, 0x50f, 0x605, 0x70c,
  66   0x80c, 0x905, 0xa0f, 0xb06, 0xc0a, 0xd03, 0xe09, 0xf00,
  67   0x190, 0x99 , 0x393, 0x29a, 0x596, 0x49f, 0x795, 0x69c,
  68   0x99c, 0x895, 0xb9f, 0xa96, 0xd9a, 0xc93, 0xf99, 0xe90,
  69   0x230, 0x339, 0x33 , 0x13a, 0x636, 0x73f, 0x435, 0x53c,
  70   0xa3c, 0xb35, 0x83f, 0x936, 0xe3a, 0xf33, 0xc39, 0xd30,
  71   0x3a0, 0x2a9, 0x1a3, 0xaa , 0x7a6, 0x6af, 0x5a5, 0x4ac,
  72   0xbac, 0xaa5, 0x9af, 0x8a6, 0xfaa, 0xea3, 0xda9, 0xca0,
  73   0x460, 0x569, 0x663, 0x76a, 0x66 , 0x16f, 0x265, 0x36c,
  74   0xc6c, 0xd65, 0xe6f, 0xf66, 0x86a, 0x963, 0xa69, 0xb60,
  75   0x5f0, 0x4f9, 0x7f3, 0x6fa, 0x1f6, 0xff , 0x3f5, 0x2fc,
  76   0xdfc, 0xcf5, 0xfff, 0xef6, 0x9fa, 0x8f3, 0xbf9, 0xaf0,
  77   0x650, 0x759, 0x453, 0x55a, 0x256, 0x35f, 0x55 , 0x15c,
  78   0xe5c, 0xf55, 0xc5f, 0xd56, 0xa5a, 0xb53, 0x859, 0x950,
  79   0x7c0, 0x6c9, 0x5c3, 0x4ca, 0x3c6, 0x2cf, 0x1c5, 0xcc ,
  80   0xfcc, 0xec5, 0xdcf, 0xcc6, 0xbca, 0xac3, 0x9c9, 0x8c0,
  81   0x8c0, 0x9c9, 0xac3, 0xbca, 0xcc6, 0xdcf, 0xec5, 0xfcc,
  82   0xcc , 0x1c5, 0x2cf, 0x3c6, 0x4ca, 0x5c3, 0x6c9, 0x7c0,
  83   0x950, 0x859, 0xb53, 0xa5a, 0xd56, 0xc5f, 0xf55, 0xe5c,
  84   0x15c, 0x55 , 0x35f, 0x256, 0x55a, 0x453, 0x759, 0x650,
  85   0xaf0, 0xbf9, 0x8f3, 0x9fa, 0xef6, 0xfff, 0xcf5, 0xdfc,
  86   0x2fc, 0x3f5, 0xff , 0x1f6, 0x6fa, 0x7f3, 0x4f9, 0x5f0,
  87   0xb60, 0xa69, 0x963, 0x86a, 0xf66, 0xe6f, 0xd65, 0xc6c,
  88   0x36c, 0x265, 0x16f, 0x66 , 0x76a, 0x663, 0x569, 0x460,
  89   0xca0, 0xda9, 0xea3, 0xfaa, 0x8a6, 0x9af, 0xaa5, 0xbac,
  90   0x4ac, 0x5a5, 0x6af, 0x7a6, 0xaa , 0x1a3, 0x2a9, 0x3a0,
  91   0xd30, 0xc39, 0xf33, 0xe3a, 0x936, 0x83f, 0xb35, 0xa3c,
  92   0x53c, 0x435, 0x73f, 0x636, 0x13a, 0x33 , 0x339, 0x230,
  93   0xe90, 0xf99, 0xc93, 0xd9a, 0xa96, 0xb9f, 0x895, 0x99c,
  94   0x69c, 0x795, 0x49f, 0x596, 0x29a, 0x393, 0x99 , 0x190,
  95   0xf00, 0xe09, 0xd03, 0xc0a, 0xb06, 0xa0f, 0x905, 0x80c,
  96   0x70c, 0x605, 0x50f, 0x406, 0x30a, 0x203, 0x109, 0x0
  97 };
  98
  99 static const int triTable[256][16] = {
 100   {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 101   { 0,  8,  3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 102   { 0,  1,  9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 103   { 1,  8,  3,  9,  8,  1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 104   { 1,  2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 105   { 0,  8,  3,  1,  2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 106   { 9,  2, 10,  0,  2,  9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 107   { 2,  8,  3,  2, 10,  8, 10,  9,  8, -1, -1, -1, -1, -1, -1, -1},
 108   { 3, 11,  2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 109   { 0, 11,  2,  8, 11,  0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 110   { 1,  9,  0,  2,  3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 111   { 1, 11,  2,  1,  9, 11,  9,  8, 11, -1, -1, -1, -1, -1, -1, -1},
 112   { 3, 10,  1, 11, 10,  3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 113   { 0, 10,  1,  0,  8, 10,  8, 11, 10, -1, -1, -1, -1, -1, -1, -1},
 114   { 3,  9,  0,  3, 11,  9, 11, 10,  9, -1, -1, -1, -1, -1, -1, -1},
 115   { 9,  8, 10, 10,  8, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 116   { 4,  7,  8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 117   { 4,  3,  0,  7,  3,  4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 118   { 0,  1,  9,  8,  4,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 119   { 4,  1,  9,  4,  7,  1,  7,  3,  1, -1, -1, -1, -1, -1, -1, -1},
 120   { 1,  2, 10,  8,  4,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 121   { 3,  4,  7,  3,  0,  4,  1,  2, 10, -1, -1, -1, -1, -1, -1, -1},
 122   { 9,  2, 10,  9,  0,  2,  8,  4,  7, -1, -1, -1, -1, -1, -1, -1},
 123   { 2, 10,  9,  2,  9,  7,  2,  7,  3,  7,  9,  4, -1, -1, -1, -1},
 124   { 8,  4,  7,  3, 11,  2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 125   {11,  4,  7, 11,  2,  4,  2,  0,  4, -1, -1, -1, -1, -1, -1, -1},
 126   { 9,  0,  1,  8,  4,  7,  2,  3, 11, -1, -1, -1, -1, -1, -1, -1},
 127   { 4,  7, 11,  9,  4, 11,  9, 11,  2,  9,  2,  1, -1, -1, -1, -1},
 128   { 3, 10,  1,  3, 11, 10,  7,  8,  4, -1, -1, -1, -1, -1, -1, -1},
 129   { 1, 11, 10,  1,  4, 11,  1,  0,  4,  7, 11,  4, -1, -1, -1, -1},
 130   { 4,  7,  8,  9,  0, 11,  9, 11, 10, 11,  0,  3, -1, -1, -1, -1},
 131   { 4,  7, 11,  4, 11,  9,  9, 11, 10, -1, -1, -1, -1, -1, -1, -1},
 132   { 9,  5,  4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 133   { 9,  5,  4,  0,  8,  3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 134   { 0,  5,  4,  1,  5,  0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 135   { 8,  5,  4,  8,  3,  5,  3,  1,  5, -1, -1, -1, -1, -1, -1, -1},
 136   { 1,  2, 10,  9,  5,  4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 137   { 3,  0,  8,  1,  2, 10,  4,  9,  5, -1, -1, -1, -1, -1, -1, -1},
 138   { 5,  2, 10,  5,  4,  2,  4,  0,  2, -1, -1, -1, -1, -1, -1, -1},
 139   { 2, 10,  5,  3,  2,  5,  3,  5,  4,  3,  4,  8, -1, -1, -1, -1},
 140   { 9,  5,  4,  2,  3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 141   { 0, 11,  2,  0,  8, 11,  4,  9,  5, -1, -1, -1, -1, -1, -1, -1},
 142   { 0,  5,  4,  0,  1,  5,  2,  3, 11, -1, -1, -1, -1, -1, -1, -1},
 143   { 2,  1,  5,  2,  5,  8,  2,  8, 11,  4,  8,  5, -1, -1, -1, -1},
 144   {10,  3, 11, 10,  1,  3,  9,  5,  4, -1, -1, -1, -1, -1, -1, -1},
 145   { 4,  9,  5,  0,  8,  1,  8, 10,  1,  8, 11, 10, -1, -1, -1, -1},
 146   { 5,  4,  0,  5,  0, 11,  5, 11, 10, 11,  0,  3, -1, -1, -1, -1},
 147   { 5,  4,  8,  5,  8, 10, 10,  8, 11, -1, -1, -1, -1, -1, -1, -1},
 148   { 9,  7,  8,  5,  7,  9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 149   { 9,  3,  0,  9,  5,  3,  5,  7,  3, -1, -1, -1, -1, -1, -1, -1},
 150   { 0,  7,  8,  0,  1,  7,  1,  5,  7, -1, -1, -1, -1, -1, -1, -1},
 151   { 1,  5,  3,  3,  5,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 152   { 9,  7,  8,  9,  5,  7, 10,  1,  2, -1, -1, -1, -1, -1, -1, -1},
 153   {10,  1,  2,  9,  5,  0,  5,  3,  0,  5,  7,  3, -1, -1, -1, -1},
 154   { 8,  0,  2,  8,  2,  5,  8,  5,  7, 10,  5,  2, -1, -1, -1, -1},
 155   { 2, 10,  5,  2,  5,  3,  3,  5,  7, -1, -1, -1, -1, -1, -1, -1},
 156   { 7,  9,  5,  7,  8,  9,  3, 11,  2, -1, -1, -1, -1, -1, -1, -1},
 157   { 9,  5,  7,  9,  7,  2,  9,  2,  0,  2,  7, 11, -1, -1, -1, -1},
 158   { 2,  3, 11,  0,  1,  8,  1,  7,  8,  1,  5,  7, -1, -1, -1, -1},
 159   {11,  2,  1, 11,  1,  7,  7,  1,  5, -1, -1, -1, -1, -1, -1, -1},
 160   { 9,  5,  8,  8,  5,  7, 10,  1,  3, 10,  3, 11, -1, -1, -1, -1},
 161   { 5,  7,  0,  5,  0,  9,  7, 11,  0,  1,  0, 10, 11, 10,  0, -1},
 162   {11, 10,  0, 11,  0,  3, 10,  5,  0,  8,  0,  7,  5,  7,  0, -1},
 163   {11, 10,  5,  7, 11,  5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 164   {10,  6,  5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 165   { 0,  8,  3,  5, 10,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 166   { 9,  0,  1,  5, 10,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 167   { 1,  8,  3,  1,  9,  8,  5, 10,  6, -1, -1, -1, -1, -1, -1, -1},
 168   { 1,  6,  5,  2,  6,  1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 169   { 1,  6,  5,  1,  2,  6,  3,  0,  8, -1, -1, -1, -1, -1, -1, -1},
 170   { 9,  6,  5,  9,  0,  6,  0,  2,  6, -1, -1, -1, -1, -1, -1, -1},
 171   { 5,  9,  8,  5,  8,  2,  5,  2,  6,  3,  2,  8, -1, -1, -1, -1},
 172   { 2,  3, 11, 10,  6,  5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 173   {11,  0,  8, 11,  2,  0, 10,  6,  5, -1, -1, -1, -1, -1, -1, -1},
 174   { 0,  1,  9,  2,  3, 11,  5, 10,  6, -1, -1, -1, -1, -1, -1, -1},
 175   { 5, 10,  6,  1,  9,  2,  9, 11,  2,  9,  8, 11, -1, -1, -1, -1},
 176   { 6,  3, 11,  6,  5,  3,  5,  1,  3, -1, -1, -1, -1, -1, -1, -1},
 177   { 0,  8, 11,  0, 11,  5,  0,  5,  1,  5, 11,  6, -1, -1, -1, -1},
 178   { 3, 11,  6,  0,  3,  6,  0,  6,  5,  0,  5,  9, -1, -1, -1, -1},
 179   { 6,  5,  9,  6,  9, 11, 11,  9,  8, -1, -1, -1, -1, -1, -1, -1},
 180   { 5, 10,  6,  4,  7,  8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 181   { 4,  3,  0,  4,  7,  3,  6,  5, 10, -1, -1, -1, -1, -1, -1, -1},
 182   { 1,  9,  0,  5, 10,  6,  8,  4,  7, -1, -1, -1, -1, -1, -1, -1},
 183   {10,  6,  5,  1,  9,  7,  1,  7,  3,  7,  9,  4, -1, -1, -1, -1},
 184   { 6,  1,  2,  6,  5,  1,  4,  7,  8, -1, -1, -1, -1, -1, -1, -1},
 185   { 1,  2,  5,  5,  2,  6,  3,  0,  4,  3,  4,  7, -1, -1, -1, -1},
 186   { 8,  4,  7,  9,  0,  5,  0,  6,  5,  0,  2,  6, -1, -1, -1, -1},
 187   { 7,  3,  9,  7,  9,  4,  3,  2,  9,  5,  9,  6,  2,  6,  9, -1},
 188   { 3, 11,  2,  7,  8,  4, 10,  6,  5, -1, -1, -1, -1, -1, -1, -1},
 189   { 5, 10,  6,  4,  7,  2,  4,  2,  0,  2,  7, 11, -1, -1, -1, -1},
 190   { 0,  1,  9,  4,  7,  8,  2,  3, 11,  5, 10,  6, -1, -1, -1, -1},
 191   { 9,  2,  1,  9, 11,  2,  9,  4, 11,  7, 11,  4,  5, 10,  6, -1},
 192   { 8,  4,  7,  3, 11,  5,  3,  5,  1,  5, 11,  6, -1, -1, -1, -1},
 193   { 5,  1, 11,  5, 11,  6,  1,  0, 11,  7, 11,  4,  0,  4, 11, -1},
 194   { 0,  5,  9,  0,  6,  5,  0,  3,  6, 11,  6,  3,  8,  4,  7, -1},
 195   { 6,  5,  9,  6,  9, 11,  4,  7,  9,  7, 11,  9, -1, -1, -1, -1},
 196   {10,  4,  9,  6,  4, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 197   { 4, 10,  6,  4,  9, 10,  0,  8,  3, -1, -1, -1, -1, -1, -1, -1},
 198   {10,  0,  1, 10,  6,  0,  6,  4,  0, -1, -1, -1, -1, -1, -1, -1},
 199   { 8,  3,  1,  8,  1,  6,  8,  6,  4,  6,  1, 10, -1, -1, -1, -1},
 200   { 1,  4,  9,  1,  2,  4,  2,  6,  4, -1, -1, -1, -1, -1, -1, -1},
 201   { 3,  0,  8,  1,  2,  9,  2,  4,  9,  2,  6,  4, -1, -1, -1, -1},
 202   { 0,  2,  4,  4,  2,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 203   { 8,  3,  2,  8,  2,  4,  4,  2,  6, -1, -1, -1, -1, -1, -1, -1},
 204   {10,  4,  9, 10,  6,  4, 11,  2,  3, -1, -1, -1, -1, -1, -1, -1},
 205   { 0,  8,  2,  2,  8, 11,  4,  9, 10,  4, 10,  6, -1, -1, -1, -1},
 206   { 3, 11,  2,  0,  1,  6,  0,  6,  4,  6,  1, 10, -1, -1, -1, -1},
 207   { 6,  4,  1,  6,  1, 10,  4,  8,  1,  2,  1, 11,  8, 11,  1, -1},
 208   { 9,  6,  4,  9,  3,  6,  9,  1,  3, 11,  6,  3, -1, -1, -1, -1},
 209   { 8, 11,  1,  8,  1,  0, 11,  6,  1,  9,  1,  4,  6,  4,  1, -1},
 210   { 3, 11,  6,  3,  6,  0,  0,  6,  4, -1, -1, -1, -1, -1, -1, -1},
 211   { 6,  4,  8, 11,  6,  8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 212   { 7, 10,  6,  7,  8, 10,  8,  9, 10, -1, -1, -1, -1, -1, -1, -1},
 213   { 0,  7,  3,  0, 10,  7,  0,  9, 10,  6,  7, 10, -1, -1, -1, -1},
 214   {10,  6,  7,  1, 10,  7,  1,  7,  8,  1,  8,  0, -1, -1, -1, -1},
 215   {10,  6,  7, 10,  7,  1,  1,  7,  3, -1, -1, -1, -1, -1, -1, -1},
 216   { 1,  2,  6,  1,  6,  8,  1,  8,  9,  8,  6,  7, -1, -1, -1, -1},
 217   { 2,  6,  9,  2,  9,  1,  6,  7,  9,  0,  9,  3,  7,  3,  9, -1},
 218   { 7,  8,  0,  7,  0,  6,  6,  0,  2, -1, -1, -1, -1, -1, -1, -1},
 219   { 7,  3,  2,  6,  7,  2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 220   { 2,  3, 11, 10,  6,  8, 10,  8,  9,  8,  6,  7, -1, -1, -1, -1},
 221   { 2,  0,  7,  2,  7, 11,  0,  9,  7,  6,  7, 10,  9, 10,  7, -1},
 222   { 1,  8,  0,  1,  7,  8,  1, 10,  7,  6,  7, 10,  2,  3, 11, -1},
 223   {11,  2,  1, 11,  1,  7, 10,  6,  1,  6,  7,  1, -1, -1, -1, -1},
 224   { 8,  9,  6,  8,  6,  7,  9,  1,  6, 11,  6,  3,  1,  3,  6, -1},
 225   { 0,  9,  1, 11,  6,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 226   { 7,  8,  0,  7,  0,  6,  3, 11,  0, 11,  6,  0, -1, -1, -1, -1},
 227   { 7, 11,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 228   { 7,  6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 229   { 3,  0,  8, 11,  7,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 230   { 0,  1,  9, 11,  7,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 231   { 8,  1,  9,  8,  3,  1, 11,  7,  6, -1, -1, -1, -1, -1, -1, -1},
 232   {10,  1,  2,  6, 11,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 233   { 1,  2, 10,  3,  0,  8,  6, 11,  7, -1, -1, -1, -1, -1, -1, -1},
 234   { 2,  9,  0,  2, 10,  9,  6, 11,  7, -1, -1, -1, -1, -1, -1, -1},
 235   { 6, 11,  7,  2, 10,  3, 10,  8,  3, 10,  9,  8, -1, -1, -1, -1},
 236   { 7,  2,  3,  6,  2,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 237   { 7,  0,  8,  7,  6,  0,  6,  2,  0, -1, -1, -1, -1, -1, -1, -1},
 238   { 2,  7,  6,  2,  3,  7,  0,  1,  9, -1, -1, -1, -1, -1, -1, -1},
 239   { 1,  6,  2,  1,  8,  6,  1,  9,  8,  8,  7,  6, -1, -1, -1, -1},
 240   {10,  7,  6, 10,  1,  7,  1,  3,  7, -1, -1, -1, -1, -1, -1, -1},
 241   {10,  7,  6,  1,  7, 10,  1,  8,  7,  1,  0,  8, -1, -1, -1, -1},
 242   { 0,  3,  7,  0,  7, 10,  0, 10,  9,  6, 10,  7, -1, -1, -1, -1},
 243   { 7,  6, 10,  7, 10,  8,  8, 10,  9, -1, -1, -1, -1, -1, -1, -1},
 244   { 6,  8,  4, 11,  8,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 245   { 3,  6, 11,  3,  0,  6,  0,  4,  6, -1, -1, -1, -1, -1, -1, -1},
 246   { 8,  6, 11,  8,  4,  6,  9,  0,  1, -1, -1, -1, -1, -1, -1, -1},
 247   { 9,  4,  6,  9,  6,  3,  9,  3,  1, 11,  3,  6, -1, -1, -1, -1},
 248   { 6,  8,  4,  6, 11,  8,  2, 10,  1, -1, -1, -1, -1, -1, -1, -1},
 249   { 1,  2, 10,  3,  0, 11,  0,  6, 11,  0,  4,  6, -1, -1, -1, -1},
 250   { 4, 11,  8,  4,  6, 11,  0,  2,  9,  2, 10,  9, -1, -1, -1, -1},
 251   {10,  9,  3, 10,  3,  2,  9,  4,  3, 11,  3,  6,  4,  6,  3, -1},
 252   { 8,  2,  3,  8,  4,  2,  4,  6,  2, -1, -1, -1, -1, -1, -1, -1},
 253   { 0,  4,  2,  4,  6,  2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 254   { 1,  9,  0,  2,  3,  4,  2,  4,  6,  4,  3,  8, -1, -1, -1, -1},
 255   { 1,  9,  4,  1,  4,  2,  2,  4,  6, -1, -1, -1, -1, -1, -1, -1},
 256   { 8,  1,  3,  8,  6,  1,  8,  4,  6,  6, 10,  1, -1, -1, -1, -1},
 257   {10,  1,  0, 10,  0,  6,  6,  0,  4, -1, -1, -1, -1, -1, -1, -1},
 258   { 4,  6,  3,  4,  3,  8,  6, 10,  3,  0,  3,  9, 10,  9,  3, -1},
 259   {10,  9,  4,  6, 10,  4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 260   { 4,  9,  5,  7,  6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 261   { 0,  8,  3,  4,  9,  5, 11,  7,  6, -1, -1, -1, -1, -1, -1, -1},
 262   { 5,  0,  1,  5,  4,  0,  7,  6, 11, -1, -1, -1, -1, -1, -1, -1},
 263   {11,  7,  6,  8,  3,  4,  3,  5,  4,  3,  1,  5, -1, -1, -1, -1},
 264   { 9,  5,  4, 10,  1,  2,  7,  6, 11, -1, -1, -1, -1, -1, -1, -1},
 265   { 6, 11,  7,  1,  2, 10,  0,  8,  3,  4,  9,  5, -1, -1, -1, -1},
 266   { 7,  6, 11,  5,  4, 10,  4,  2, 10,  4,  0,  2, -1, -1, -1, -1},
 267   { 3,  4,  8,  3,  5,  4,  3,  2,  5, 10,  5,  2, 11,  7,  6, -1},
 268   { 7,  2,  3,  7,  6,  2,  5,  4,  9, -1, -1, -1, -1, -1, -1, -1},
 269   { 9,  5,  4,  0,  8,  6,  0,  6,  2,  6,  8,  7, -1, -1, -1, -1},
 270   { 3,  6,  2,  3,  7,  6,  1,  5,  0,  5,  4,  0, -1, -1, -1, -1},
 271   { 6,  2,  8,  6,  8,  7,  2,  1,  8,  4,  8,  5,  1,  5,  8, -1},
 272   { 9,  5,  4, 10,  1,  6,  1,  7,  6,  1,  3,  7, -1, -1, -1, -1},
 273   { 1,  6, 10,  1,  7,  6,  1,  0,  7,  8,  7,  0,  9,  5,  4, -1},
 274   { 4,  0, 10,  4, 10,  5,  0,  3, 10,  6, 10,  7,  3,  7, 10, -1},
 275   { 7,  6, 10,  7, 10,  8,  5,  4, 10,  4,  8, 10, -1, -1, -1, -1},
 276   { 6,  9,  5,  6, 11,  9, 11,  8,  9, -1, -1, -1, -1, -1, -1, -1},
 277   { 3,  6, 11,  0,  6,  3,  0,  5,  6,  0,  9,  5, -1, -1, -1, -1},
 278   { 0, 11,  8,  0,  5, 11,  0,  1,  5,  5,  6, 11, -1, -1, -1, -1},
 279   { 6, 11,  3,  6,  3,  5,  5,  3,  1, -1, -1, -1, -1, -1, -1, -1},
 280   { 1,  2, 10,  9,  5, 11,  9, 11,  8, 11,  5,  6, -1, -1, -1, -1},
 281   { 0, 11,  3,  0,  6, 11,  0,  9,  6,  5,  6,  9,  1,  2, 10, -1},
 282   {11,  8,  5, 11,  5,  6,  8,  0,  5, 10,  5,  2,  0,  2,  5, -1},
 283   { 6, 11,  3,  6,  3,  5,  2, 10,  3, 10,  5,  3, -1, -1, -1, -1},
 284   { 5,  8,  9,  5,  2,  8,  5,  6,  2,  3,  8,  2, -1, -1, -1, -1},
 285   { 9,  5,  6,  9,  6,  0,  0,  6,  2, -1, -1, -1, -1, -1, -1, -1},
 286   { 1,  5,  8,  1,  8,  0,  5,  6,  8,  3,  8,  2,  6,  2,  8, -1},
 287   { 1,  5,  6,  2,  1,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 288   { 1,  3,  6,  1,  6, 10,  3,  8,  6,  5,  6,  9,  8,  9,  6, -1},
 289   {10,  1,  0, 10,  0,  6,  9,  5,  0,  5,  6,  0, -1, -1, -1, -1},
 290   { 0,  3,  8,  5,  6, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 291   {10,  5,  6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 292   {11,  5, 10,  7,  5, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 293   {11,  5, 10, 11,  7,  5,  8,  3,  0, -1, -1, -1, -1, -1, -1, -1},
 294   { 5, 11,  7,  5, 10, 11,  1,  9,  0, -1, -1, -1, -1, -1, -1, -1},
 295   {10,  7,  5, 10, 11,  7,  9,  8,  1,  8,  3,  1, -1, -1, -1, -1},
 296   {11,  1,  2, 11,  7,  1,  7,  5,  1, -1, -1, -1, -1, -1, -1, -1},
 297   { 0,  8,  3,  1,  2,  7,  1,  7,  5,  7,  2, 11, -1, -1, -1, -1},
 298   { 9,  7,  5,  9,  2,  7,  9,  0,  2,  2, 11,  7, -1, -1, -1, -1},
 299   { 7,  5,  2,  7,  2, 11,  5,  9,  2,  3,  2,  8,  9,  8,  2, -1},
 300   { 2,  5, 10,  2,  3,  5,  3,  7,  5, -1, -1, -1, -1, -1, -1, -1},
 301   { 8,  2,  0,  8,  5,  2,  8,  7,  5, 10,  2,  5, -1, -1, -1, -1},
 302   { 9,  0,  1,  5, 10,  3,  5,  3,  7,  3, 10,  2, -1, -1, -1, -1},
 303   { 9,  8,  2,  9,  2,  1,  8,  7,  2, 10,  2,  5,  7,  5,  2, -1},
 304   { 1,  3,  5,  3,  7,  5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 305   { 0,  8,  7,  0,  7,  1,  1,  7,  5, -1, -1, -1, -1, -1, -1, -1},
 306   { 9,  0,  3,  9,  3,  5,  5,  3,  7, -1, -1, -1, -1, -1, -1, -1},
 307   { 9,  8,  7,  5,  9,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 308   { 5,  8,  4,  5, 10,  8, 10, 11,  8, -1, -1, -1, -1, -1, -1, -1},
 309   { 5,  0,  4,  5, 11,  0,  5, 10, 11, 11,  3,  0, -1, -1, -1, -1},
 310   { 0,  1,  9,  8,  4, 10,  8, 10, 11, 10,  4,  5, -1, -1, -1, -1},
 311   {10, 11,  4, 10,  4,  5, 11,  3,  4,  9,  4,  1,  3,  1,  4, -1},
 312   { 2,  5,  1,  2,  8,  5,  2, 11,  8,  4,  5,  8, -1, -1, -1, -1},
 313   { 0,  4, 11,  0, 11,  3,  4,  5, 11,  2, 11,  1,  5,  1, 11, -1},
 314   { 0,  2,  5,  0,  5,  9,  2, 11,  5,  4,  5,  8, 11,  8,  5, -1},
 315   { 9,  4,  5,  2, 11,  3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 316   { 2,  5, 10,  3,  5,  2,  3,  4,  5,  3,  8,  4, -1, -1, -1, -1},
 317   { 5, 10,  2,  5,  2,  4,  4,  2,  0, -1, -1, -1, -1, -1, -1, -1},
 318   { 3, 10,  2,  3,  5, 10,  3,  8,  5,  4,  5,  8,  0,  1,  9, -1},
 319   { 5, 10,  2,  5,  2,  4,  1,  9,  2,  9,  4,  2, -1, -1, -1, -1},
 320   { 8,  4,  5,  8,  5,  3,  3,  5,  1, -1, -1, -1, -1, -1, -1, -1},
 321   { 0,  4,  5,  1,  0,  5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 322   { 8,  4,  5,  8,  5,  3,  9,  0,  5,  0,  3,  5, -1, -1, -1, -1},
 323   { 9,  4,  5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 324   { 4, 11,  7,  4,  9, 11,  9, 10, 11, -1, -1, -1, -1, -1, -1, -1},
 325   { 0,  8,  3,  4,  9,  7,  9, 11,  7,  9, 10, 11, -1, -1, -1, -1},
 326   { 1, 10, 11,  1, 11,  4,  1,  4,  0,  7,  4, 11, -1, -1, -1, -1},
 327   { 3,  1,  4,  3,  4,  8,  1, 10,  4,  7,  4, 11, 10, 11,  4, -1},
 328   { 4, 11,  7,  9, 11,  4,  9,  2, 11,  9,  1,  2, -1, -1, -1, -1},
 329   { 9,  7,  4,  9, 11,  7,  9,  1, 11,  2, 11,  1,  0,  8,  3, -1},
 330   {11,  7,  4, 11,  4,  2,  2,  4,  0, -1, -1, -1, -1, -1, -1, -1},
 331   {11,  7,  4, 11,  4,  2,  8,  3,  4,  3,  2,  4, -1, -1, -1, -1},
 332   { 2,  9, 10,  2,  7,  9,  2,  3,  7,  7,  4,  9, -1, -1, -1, -1},
 333   { 9, 10,  7,  9,  7,  4, 10,  2,  7,  8,  7,  0,  2,  0,  7, -1},
 334   { 3,  7, 10,  3, 10,  2,  7,  4, 10,  1, 10,  0,  4,  0, 10, -1},
 335   { 1, 10,  2,  8,  7,  4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 336   { 4,  9,  1,  4,  1,  7,  7,  1,  3, -1, -1, -1, -1, -1, -1, -1},
 337   { 4,  9,  1,  4,  1,  7,  0,  8,  1,  8,  7,  1, -1, -1, -1, -1},
 338   { 4,  0,  3,  7,  4,  3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 339   { 4,  8,  7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 340   { 9, 10,  8, 10, 11,  8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 341   { 3,  0,  9,  3,  9, 11, 11,  9, 10, -1, -1, -1, -1, -1, -1, -1},
 342   { 0,  1, 10,  0, 10,  8,  8, 10, 11, -1, -1, -1, -1, -1, -1, -1},
 343   { 3,  1, 10, 11,  3, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 344   { 1,  2, 11,  1, 11,  9,  9, 11,  8, -1, -1, -1, -1, -1, -1, -1},
 345   { 3,  0,  9,  3,  9, 11,  1,  2,  9,  2, 11,  9, -1, -1, -1, -1},
 346   { 0,  2, 11,  8,  0, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 347   { 3,  2, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 348   { 2,  3,  8,  2,  8, 10, 10,  8,  9, -1, -1, -1, -1, -1, -1, -1},
 349   { 9, 10,  2,  0,  9,  2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 350   { 2,  3,  8,  2,  8, 10,  0,  1,  8,  1, 10,  8, -1, -1, -1, -1},
 351   { 1, 10,  2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 352   { 1,  3,  8,  9,  1,  8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 353   { 0,  9,  1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 354   { 0,  3,  8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
 355   {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}
 356 };
 357
 358
 359
 360 /* Linearly interpolate the position where an isosurface cuts
 361    an edge between two vertices, each with their own scalar value
 362 */
 363 static XYZ
 364 interp_vertex (double isolevel, XYZ p1, XYZ p2, double valp1, double valp2)
 365 {
 366   double mu;
 367   XYZ p;
 368
 369   if (ABS(isolevel-valp1) < 0.00001)
 370     return(p1);
 371   if (ABS(isolevel-valp2) < 0.00001)
 372     return(p2);
 373   if (ABS(valp1-valp2) < 0.00001)
 374     return(p1);
 375   mu = (isolevel - valp1) / (valp2 - valp1);
 376   p.x = p1.x + mu * (p2.x - p1.x);
 377   p.y = p1.y + mu * (p2.y - p1.y);
 378   p.z = p1.z + mu * (p2.z - p1.z);
 379
 380   return(p);
 381 }
 382
 383
 384 /* Given a grid cell and an isolevel, calculate the triangular
 385    facets required to represent the isosurface through the cell.
 386    Return the number of triangular facets.
 387    `triangles' will be loaded up with the vertices at most 5 triangular facets.
 388    0 will be returned if the grid cell is either totally above
 389    of totally below the isolevel.
 390
 391   By Paul Bourke <pbourke@swin.edu.au>
 392 */
 393 static int
 394 march_one_cube (GRIDCELL grid, double isolevel, TRIANGLE *triangles)
 395 {
 396   int i, ntriang;
 397   int cubeindex;
 398   XYZ vertlist[12];
 399
 400   /*
 401     Determine the index into the edge table which
 402     tells us which vertices are inside of the surface
 403   */
 404   cubeindex = 0;
 405   if (grid.val[0] < isolevel) cubeindex |= 1;
 406   if (grid.val[1] < isolevel) cubeindex |= 2;
 407   if (grid.val[2] < isolevel) cubeindex |= 4;
 408   if (grid.val[3] < isolevel) cubeindex |= 8;
 409   if (grid.val[4] < isolevel) cubeindex |= 16;
 410   if (grid.val[5] < isolevel) cubeindex |= 32;
 411   if (grid.val[6] < isolevel) cubeindex |= 64;
 412   if (grid.val[7] < isolevel) cubeindex |= 128;
 413
 414   /* Cube is entirely in/out of the surface */
 415   if (edgeTable[cubeindex] == 0)
 416     return(0);
 417
 418   /* Find the vertices where the surface intersects the cube */
 419   if (edgeTable[cubeindex] & 1)
 420     vertlist[0] =
 421       interp_vertex (isolevel,grid.p[0],grid.p[1],grid.val[0],grid.val[1]);
 422   if (edgeTable[cubeindex] & 2)
 423     vertlist[1] =
 424       interp_vertex (isolevel,grid.p[1],grid.p[2],grid.val[1],grid.val[2]);
 425   if (edgeTable[cubeindex] & 4)
 426     vertlist[2] =
 427       interp_vertex (isolevel,grid.p[2],grid.p[3],grid.val[2],grid.val[3]);
 428   if (edgeTable[cubeindex] & 8)
 429     vertlist[3] =
 430       interp_vertex (isolevel,grid.p[3],grid.p[0],grid.val[3],grid.val[0]);
 431   if (edgeTable[cubeindex] & 16)
 432     vertlist[4] =
 433       interp_vertex (isolevel,grid.p[4],grid.p[5],grid.val[4],grid.val[5]);
 434   if (edgeTable[cubeindex] & 32)
 435     vertlist[5] =
 436       interp_vertex (isolevel,grid.p[5],grid.p[6],grid.val[5],grid.val[6]);
 437   if (edgeTable[cubeindex] & 64)
 438     vertlist[6] =
 439       interp_vertex (isolevel,grid.p[6],grid.p[7],grid.val[6],grid.val[7]);
 440   if (edgeTable[cubeindex] & 128)
 441     vertlist[7] =
 442       interp_vertex (isolevel,grid.p[7],grid.p[4],grid.val[7],grid.val[4]);
 443   if (edgeTable[cubeindex] & 256)
 444     vertlist[8] =
 445       interp_vertex (isolevel,grid.p[0],grid.p[4],grid.val[0],grid.val[4]);
 446   if (edgeTable[cubeindex] & 512)
 447     vertlist[9] =
 448       interp_vertex (isolevel,grid.p[1],grid.p[5],grid.val[1],grid.val[5]);
 449   if (edgeTable[cubeindex] & 1024)
 450     vertlist[10] =
 451       interp_vertex (isolevel,grid.p[2],grid.p[6],grid.val[2],grid.val[6]);
 452   if (edgeTable[cubeindex] & 2048)
 453     vertlist[11] =
 454       interp_vertex (isolevel,grid.p[3],grid.p[7],grid.val[3],grid.val[7]);
 455
 456   /* Create the triangle */
 457   ntriang = 0;
 458   for (i=0; triTable[cubeindex][i] != -1; i+=3)
 459     {
 460       triangles[ntriang].p[0] = vertlist[triTable[cubeindex][i  ]];
 461       triangles[ntriang].p[1] = vertlist[triTable[cubeindex][i+1]];
 462       triangles[ntriang].p[2] = vertlist[triTable[cubeindex][i+2]];
 463       ntriang++;
 464     }
 465
 466   return(ntriang);
 467 }
 468
 469
 470 /* Walking the grid.  By jwz.
 471  */
 472
 473
 474 /* Computes the normal of the scalar field at the given point,
 475    for vertex normals (as opposed to face normals.)
 476  */
 477 static void
 478 do_function_normal (double x, double y, double z,
 479                     double (*compute_fn) (double x, double y, double z,
 480                                           void *closure),
 481                     void *c)
 482 {
 483   XYZ n;
 484   double off = 0.5;
 485   n.x = compute_fn (x-off, y, z, c) - compute_fn (x+off, y, z, c);
 486   n.y = compute_fn (x, y-off, z, c) - compute_fn (x, y+off, z, c);
 487   n.z = compute_fn (x, y, z-off, c) - compute_fn (x, y, z+off, c);
 488   /* normalize (&n); */
 489   glNormal3f (n.x, n.y, n.z);
 490 }
 491
 492
 493 /* Given a function capable of generating a value at any XYZ position,
 494    creates OpenGL faces for the solids defined.
 495
 496    init_fn is called at the beginning for initial, and returns an object.
 497    free_fn is called at the end.
 498
 499    compute_fn is called for each XYZ in the specified grid, and returns
 500    the double value of that coordinate.  If smoothing is on, then
 501    compute_fn will also be called twice more for each emitted vertex,
 502    in order to calculate vertex normals (so don't assume it will only
 503    be called with values falling on the grid boundaries.)
 504
 505    Points are inside an object if the are less than `isolevel', and
 506    outside otherwise.
 507 */
 508 void
 509 marching_cubes (int grid_size,     /* density of the mesh */
 510                 double isolevel,   /* cutoff point for "in" versus "out" */
 511                 int wireframe_p,   /* wireframe, or solid */
 512                 int smooth_p,      /* smooth, or faceted */
 513
 514                 void * (*init_fn)    (double grid_size, void *closure1),
 515                 double (*compute_fn) (double x, double y, double z,
 516                                       void *closure2),
 517                 void   (*free_fn)    (void *closure2),
 518                 void *closure1,
 519
 520                 unsigned long *polygon_count)
 521 {
 522   int planesize = grid_size * grid_size;
 523   int x, y, z;
 524   void *closure2 = 0;
 525   unsigned long polys = 0;
 526   double *layers;
 527
 528   layers = (double *) calloc (sizeof (*layers), planesize * 2);
 529   if (!layers)
 530     {
 531       fprintf (stderr, "%s: out of memory for %dx%dx%d grid\n",
 532                progname, grid_size, grid_size, 2);
 533       exit (1);
 534     }
 535
 536   if (init_fn)
 537     closure2 = init_fn (grid_size, closure1);
 538
 539   glFrontFace(GL_CCW);
 540   if (!wireframe_p)
 541     glBegin (GL_TRIANGLES);
 542
 543   for (z = 0; z < grid_size; z++)
 544     {
 545       double *layer0 = (z & 1 ? layers+planesize : layers);
 546       double *layer1 = (z & 1 ? layers           : layers+planesize);
 547       double *row;
 548
 549       /* Fill in the XY grid on the currently-bottommost layer. */
 550       row = layer1;
 551       for (y = 0; y < grid_size; y++, row += grid_size)
 552         {
 553           double *cell = row;
 554           for (x = 0; x < grid_size; x++, cell++)
 555             *cell = compute_fn (x, y, z, closure2);
 556         }
 557
 558       /* Now we've completed one layer (an XY slice of Z.)  Now we can
 559          generate the polygons that fill the space between this layer
 560          and the previous one (unless this is the first layer.)
 561       */
 562       if (z == 0) continue;
 563
 564       for (y = 1; y < grid_size; y += 1)
 565         for (x = 1; x < grid_size; x += 1)
 566           {
 567             TRIANGLE tri[6];
 568             int i, ntri;
 569             GRIDCELL cell;
 570
 571             /* This is kinda hokey, there ought to be a more efficient
 572                way to do this... */
 573             cell.p[0].x = x-1; cell.p[0].y = y-1; cell.p[0].z = z-1;
 574             cell.p[1].x = x  ; cell.p[1].y = y-1; cell.p[1].z = z-1;
 575             cell.p[2].x = x  ; cell.p[2].y = y  ; cell.p[2].z = z-1;
 576             cell.p[3].x = x-1; cell.p[3].y = y  ; cell.p[3].z = z-1;
 577             cell.p[4].x = x-1; cell.p[4].y = y-1; cell.p[4].z = z  ;
 578             cell.p[5].x = x  ; cell.p[5].y = y-1; cell.p[5].z = z  ;
 579             cell.p[6].x = x  ; cell.p[6].y = y  ; cell.p[6].z = z  ;
 580             cell.p[7].x = x-1; cell.p[7].y = y  ; cell.p[7].z = z  ;
 581
 582 # define GRID(X,Y,WHICH) ((WHICH) \
 583                           ? layer1[((Y)*grid_size) + ((X))] \
 584                           : layer0[((Y)*grid_size) + ((X))])
 585
 586             cell.val[0] = GRID (x-1, y-1, 0);
 587             cell.val[1] = GRID (x  , y-1, 0);
 588             cell.val[2] = GRID (x  , y  , 0);
 589             cell.val[3] = GRID (x-1, y  , 0);
 590             cell.val[4] = GRID (x-1, y-1, 1);
 591             cell.val[5] = GRID (x  , y-1, 1);
 592             cell.val[6] = GRID (x  , y  , 1);
 593             cell.val[7] = GRID (x-1, y  , 1);
 594 # undef GRID
 595
 596             /* Now generate the triangles for this cubic segment,
 597                and emit the GL faces.
 598             */
 599             ntri = march_one_cube (cell, isolevel, tri);
 600             polys += ntri;
 601             for (i = 0; i < ntri; i++)
 602               {
 603                 if (wireframe_p) glBegin (GL_LINE_LOOP);
 604
 605                 /* If we're smoothing, we need to call the field function
 606                    again for each vertex (via function_normal().)  If we're
 607                    not smoothing, then we can just compute the normal from
 608                    this triangle.
 609                  */
 610                 if (!smooth_p)
 611                   do_normal (tri[i].p[0].x, tri[i].p[0].y, tri[i].p[0].z,
 612                              tri[i].p[1].x, tri[i].p[1].y, tri[i].p[1].z,
 613                              tri[i].p[2].x, tri[i].p[2].y, tri[i].p[2].z);
 614
 615 # define VERT(X,Y,Z) \
 616                 if (smooth_p) \
 617                   do_function_normal ((X), (Y), (Z), compute_fn, closure2); \
 618                 glVertex3f ((X), (Y), (Z))
 619
 620                 VERT (tri[i].p[0].x, tri[i].p[0].y, tri[i].p[0].z);
 621                 VERT (tri[i].p[1].x, tri[i].p[1].y, tri[i].p[1].z);
 622                 VERT (tri[i].p[2].x, tri[i].p[2].y, tri[i].p[2].z);
 623 # undef VERT
 624                 if (wireframe_p) glEnd ();
 625               }
 626           }
 627     }
 628
 629   if (!wireframe_p)
 630     glEnd ();
 631
 632   free (layers);
 633
 634   if (free_fn)
 635     free_fn (closure2);
 636
 637   if (polygon_count)
 638     *polygon_count = polys;
 639 }