GRAPP 2007
Adaptive Cube Tessellation for Topologically Correct Isosurfaces
First presented at the International Conference on Computer Graphics Theory and Applications (GRAPP) 2007,
extended and revised for JVRB
urn:nbn:de:0009613098
Abstract
Three dimensional datasets representing scalar fields are frequently rendered using isosurfaces. For datasets arranged as a cubic lattice, the marching cubes algorithm is the most used isosurface extraction method. However, the marching cubes algorithm produces some ambiguities which have been solved using different approaches that normally imply a more complex process. One of them is to tessellate the cubes into tetrahedra, and by using a similar method (marching tetrahedra), to build the isosurface. The main drawback of other tessellations is that they do not produce the same isosurface topologies as those generated by improved marching cubes algorithms. We propose an adaptive tessellation that, being independent of the isovalue, preserves the topology. Moreover the tessellationallows the isosurface to evolve continuously when the isovalue is changed continuously.
Keywords: Volume Visualization, Isosurfaces, Marching Cubes, Marching Tetrahedra
Subjects: Computer Graphics, Visualization, Tesselation
A popular representation of digitized volumes is a regular lattice of points which represent the digitized samples. This set of samples can be used directly or can be previously filtered to obtain a reduced dataset. The most used lattice for some fields is a grid of cubic cells. This representation is easily obtained from a computerized axial tomography (CAT scan) or by magnetic resonance imaging (MRI). In both cases several 2D images are obtained and can be easily represented by a 3D grid in which every vertex represents a point sample.
Thus, a volume can be represented by the equation:
({(x,y,z,γ)}, F(x,y,z) :
³ → Γ) (1)
where{(x,y,z,γ)} is the set of samples in the property domain Γ, and F(x,y,z) : ³ → Γ is the function that approximates values between samples.
There are several methods to visualize a volume represented by a grid: by slices, direct visualization [ Lev90 ], or by means of isosurface extraction. By this later method, the focus of this paper, a threshold valueγ_{0} is set and an isosurface F(x,y,z) = γ_{0} is built and rendered as a volume visualization. Usually γ_{0} is also called isovalue. Several parts of the volume can be rendered by modifying γ_{0} . In this method F(x,y,z) must be extended outside the samples, usually as an interpolation function. The classical algorithm in this category is the marching cubes method [ LC87, NY06 ]. This method builds the isosurface cube by cube, marching through all the cubes within the grid of cubic cells and extends F as the linear interpolation along the edges of every cube.
Every cell is classified, by comparing its eight vertex values with the threshold γ_{0} , as belonging to one of the fifteen possible equivalence classes. Every class (or case) has an isosurface represented as a triangle mesh. The entire isosurface to be rendered is obtained by joining the pieces of isosurface generated for every cube. Figure 1 shows the fifteen cases with their interior triangle meshes.
However, this method is ambiguous in some topologycal aspects. The ambiguity has been well studied in [ WvG90, NH91, MSS94, Che95, Lop99, CGMS00, LB03, Nie03 ]. For instance, figure 2a shows that a crack arises on the isosurface between the two cells; figure 2b shows the triangle mesh proposed by marching cubes for the Case 4 however, the one in figure 2c may also be possible. To find out which one is correct we need information about the interior of the cell.
Lopes, Brodlie and Nielson in [ LB03, Nie03 ] solve the ambiguity studying the interior of the faces and the interior of the cell using the trilinear interpolation equation to define the property variation. As a consequence, the number of different equivalence classes is increased to 31 and 57 cases respectively. In this paper we will refer to these method as extended marching cubes(EMC).
Other authors [ ZCT95, GH95, CP98, GP00 ] solve the ambiguity using a tessellation of the cell into tetrahedra and building the isosurface by marching tetrahedra [ PT90 ]. This method is similar to marching cubes but based on a tetrahedron instead of a cube, with the added advantage that, in this method, there are no ambiguities as it is well known and only has 3 equivalence classes.
However, for many classes of EMC [ LB03, Nie03 ], the topology of the isosurface which is built by these methods is not the same as the topology of the isosurface built by marching tetrahedra from tessellations in other published works.
In this paper we propose a tessellation of cells into tetrahedra which produces the same isosurface topology as the one that would be extracted by the extended marching cubes methods.
The next section analyses the sources of ambiguity in a cube. Section three presents some tessellations used to obtain tetrahedra from cubes. In section four we put forth our proposal of tessellation. The paper ends with a section presenting the results and the conclusions.
In order to determine the correct topology for an isosurface we have to study the interpolation function F, which is a trilinear interpolation function. For a more exhaustive study the paper by Nielson [ Nie03 ] may be consulted.
Let us begin analyzing a cube face. Taking into account that every face vertex can be positive (its value is greater than the isovalue) or negative and that a face has four vertices, there are 2^{4} = 16 possible configurations of a face, however by rotation or by complementation (positive/negative), there are only 4 equivalence classes (see figure 3).
Assuming a linear variation along the edges, it is easy to prove that Case 0 does not have an isocurve and Cases 1 and 2 have the two isocurve topologies shown in figure 4. However, Case 3 has two possible isocurve topologies as is shown in figure 4 depending on which vertices can be connected by a line without intersecting the isocurves, the negative ones (case 3a) or the positive ones (case 3b).
To discover which one is the correct topology we have
to compute the face saddle point, the point where one
topology changes into another. This point is computed by
using the function that interpolates the interior of the face,
which is the bilinear interpolation:
F_{face}(x,y) = axy + bx + cy + d (2)
The face saddle point is calculated by making the
partial derivates, with respect to x and y, equal to zero:
and its value is computed using the equation 2.
So, we have a 5^{th} point which can be positive or negative and allows us to choose the correct topology of the isocurve as it is shown in figure 5. The face saddle point becomes the inflection point between the two configurations; it is a contact point when the isovalue is equal to the face saddle point value (central image in figure 5).
Figure 5. An isocurve which is continuously moved when the isovalue changes. The labels '+', '' and '=' mean if the value of that point is greater, lesser or equal than the isovalue.
The face saddle point (FSP), whenever it exists and is inside the face, can be used to tessellate it into 4 triangles. So the correct topology is directly obtained by processing the triangles instead of the square (see figure 6).
If the face saddle point is not present or is outside the face, it can be shown that the face will be always classified as belonging to Cases 0, 1 or 2 (figure 3), which are not ambiguous.
For the interior of the cube the process is similar: the body saddle points need to be computed.
The function that interpolates the interior of the cube is
the trilinear interpolation:
F(x,y,z) = axyz + bxy + cyz + dzx + ex + fy + gz + h (4)
The body saddle point (BSP) is also obtained by making 0 the three partial derivates:
The result of this equation system gives two possible body saddle points which are used to solve the ambiguities and choose the correct topology for the isosurface. The body saddle points are inflection points between the different configurations for an ambiguous cell. In this way, we can use the body saddle points, when they are present and they are inside the cell, to tessellate the cube into tetrahedra. The correct topology of the isosurface can be obtained by marching tetrahedra as will be shown in sections 4 and 5.
All the possible topologies, numbered following Lopes' methodology in [ LB03 ], are shown in figure 7. The first number defines the case taking into account the configuration of positive and negative vertices, the second number defines different solutions for ambiguous faces, and the third one defines different solutions for the ambiguous body. The case 0, without isosurface, is not shown. Nielson presents more cases in [ Nie03 ], this is due to the use of a different equivalence relationship, but many of the cases are equivalents according to the Lopes' equivalence relationship.
These topologies can be easily obtained by marching tetrahedra with the adequate tessellation. The next section shows related works on tessellations and section 4 puts forth our proposal for the tessellation.
The marching tetrahedra method [ PT90 ] for isosurface building is similar to the marching cubes method [ LC87 ]:
This method is non ambiguous by assuming a linear interpolation along the edges. Thus, by tessellating a cell into tetrahedra (figure 8 shows the most used tessellations) the ambiguity problem could be solved as can be seen in [ PT90, GH95 ].
However, these studies do not take into account that the property variation along a diagonal of a cube face is not linear, but quadratic, so that diagonal which is taken as a tetrahedron edge could have two intersections with the isosurface and not only one. Zhou et al. [ ZCT95 ] study this question and define a new way to build isosurfaces inside tetrahedra taking into account the quadratic variation of properties along diagonal edges, but there are 59 different cases! And moreover the method does not produce topologies like that of 4.1.2 in figure 7.
Chan and Purisima [ CP98 ] define a different tessellation (see figure 9) by defining tetraedra between adjacent cells, but the method does not produce topologies like that of 13.5.1 in figure 7.
Later studies have improved this method. Treece et al. [ TPG99 ] reduce the number or triangles by clustering tetrahedra vertices which preserves the topology. The generated triangles are more regular, so the Gouraud shading is improved. Gerstner et al. [ GP00 ] achieve a multiresolution tessellation by a recursive bisection of some tetrahedron in two. They take into account several critical points to preserve or not the topology in accordance with certain criteria. Both studies start from a tessellation that does not produce all the topologies of figure 7.
Chiang et al. [ CL03 ] propose a progressive simplification of tetrahedral meshes by using a contour tree, a data structure that represents the relations between connected components of the isosurfaces embedded in a volume dataset. Since the structure used includes critical points, the simplification preserves the topology. The authors focus on an irregular grid, however their method can be also applied to regular grids.
Our proposal produces a tessellation that preserves all the topologies of figure 7 and can be used as an initial tessellation for the method of Chiang.
We propose to carry out an adaptive tessellation of every cell into tetrahedra in such a way that:

The tessellation must be isovalue independent in order to compute it just once, avoiding extra computations every time that the isovalue is changed.

The isosurface built inside the cell by joining the pieces of isosurfaces built from every tetrahedron of the tessellation has to be topologically equivalent to the isosurface built by using the extended marching cubes already commented in section 2 and in figure 7.

The isosurface has to move continuously when the isovalue changes continuously.
It will be done by taking into account the face and body saddle points.
We will use a basic tessellation which is valid for all non ambiguous cell configurations, that is to say, it is valid for those cells without face saddle points nor body saddle points. Then, this basic tessellation will be modified when face or body saddle points exist. These special points will be tetrahedra vertices and will be located at their exact position to improve the isosurface accuracy.
Basic case: 0 ambiguous faces and 0 body saddle points
(Case 00)
The basic case corresponds to cells without ambiguous
faces and without body saddle points. This kind of cell is
tessellated as shown on figure 8a. This tessellation
produces 6 tetrahedra. Note that it is possible to build a
tessellation which produces just 5 tetrahedra (see figure
8b) but it is less homogeneous than the one proposed,
because it implies the use of 2 possible symmetric
tessellations depending on the position of the cell as can
be deduced in figure 8b. Another advantage is that
this tesselation is similar to the others that will be
shown in the paper. In figure 10 you can see an
example of this case where the isovalue is continuously
changed.
Case 01
This case corresponds to cells without ambiguous
faces and just 1 body saddle point. This kind of cell is
tessellated as is shown in figure 11a where the body
saddle point is the central point. In order to see it clearly,
one must look at figure 11b where the tessellation
which corresponds to a face is shown. The other 5 faces
are tessellated in a similar way. The cube is tesselated
into 6 pyramids and each pyramid is tesselated into
2 tetrahedra. As such, this tessellation produces 12
tetrahedra. For clarity, in the figure 11a the body saddle
point is drawn in the center of the cell, but it will always
be positioned in its exact position. All the tetrahedra share
the body saddle point, so the topology around this point is
preserved.
Figure 11. Tessellation for the case 01. To understand the tessellation, imagine the cube composed by 6 pyramids where every pyramid has got its square base on a cube face and its top vertex on the body saddle point. Then, imagine that every pyramid is composed by 2 tetrahedra as it is shown in (b).
Cases n1 (1 ≤ n ≤ 6)
These cases corresponds to cells with just 1 body saddle
point and several ambiguous faces with one face saddle
point associated with each ambiguous face. These cases
are tesselated in a similar way to the case 01, the cube is
tesselated into 6 pyramids, one pyramid for each face. The
pyramids for the nonambiguous faces are tesselated into 2
tetrahedra, as we have shown for the case 01 (figure
11b). The pyramids for the ambiguous faces are
tesselated into 4 tetrahedra by using the face saddle point
in the base of the pyramid; as you can see in the figure
12b. The figure 12a shows the cube tesselated for the
case 11 with the ambiguous face in the bottom face
of the cube. All the saddle points will be positioned
in their exact positions. This tessellation produces
12 + 2n tetrahedra. In these cases the topology is
preserved around the face and body saddle points as can
be easily observed by the fact that all the tetrahedra of an
ambiguous face share the face saddle point and all the
tetrahedra from the tessellation share the body saddle
point.
Figure 12. Analysis for the case 11. To understand the tessellation, imagine the cube tessellated like the case 00 but for the bottom face, its pyramid is tesselated into 4 tetrahedra as it is shown in (b).
In the figure 13 you can see an example for the case 11. The isovalue changes continuously and so, the isosurface is moved continuously from one topology to another.
Figure 13. Isosurface in a cell within 1 body saddle point and 1 face saddle point it its bottom face. The right bottom image is the same case than the left bottom one, only the viewpoint changes.
Cases n0 (1 ≤ n ≤ 6)
The cases n0 are transformed into cases n1 by adding a
fictitious body saddle point in the center of the cube. We
proceed this way to have tesselations analogous to the
others. The figures 14, 15, and 16 show several
examples with different configurations.
Case 02
This represents a cell with no ambiguous faces and 2 body
saddle points. These two body saddle points form
an edge
. Using this edge and the 6 edges of
the cuboid formed by the vertices V1, V5, V4, V6,
V2, V3 (see figure 17a), 6 tetrahedra are built (see
figure 17b). The back, left and bottom faces (which
share V0) are tessellated as it was shown in figure
11b but using B1 as the pyramid vertex (the body
saddle point nearest to the cell origin). The front, right
and top faces (which share V7) are tessellated in the
same way but using B2 as the pyramid vertex, the
other body saddle point. This tessellation produces 18
tetrahedra.
Figure 17. Analysis for the case 02. To understand this figure, imagine the cube without its vertices V0 and V7. Then a tetrahedron is built among (B1, B2, V1, V5). Another one is built among (B1, B2, V5, V4). And so on with (B1, B2, V4, V6), (B1, B2, V6, V2), (B1, B2, V2, V3) and (B1, B2, V3, V1). These six tetrahedra are the figure (b).
Cases n2 (1 ≤ n ≤ 6)
The tessellation according to this case is carried out as the
one shown for the case 02 but the ambiguous faces are
tessellated as it was shown in figure 12b instead of the
tessellation shown in figure 11b. This tessellation
produces 18 + 2n tetrahedra. The topology is also
preserved in theese cases as can be seen in figure
18.
Figure 18. Isosurface in a cell with 2 body saddle points. The images in the third row show the first and last cases of the first row from other viewpoint. The image of the second row shows an intermediate case.
The main advantage of this method is that tessellations, being independent of isovalue, produce the same isosurface topologies as the ones obtained by the recents improvements of marching cubes [ LB03, Nie03 ] but with less computations because, with our method, an isovalue change does not produce a new tessellation, but rather marching tetrahedra has to be run again; whereas with Lopes and Nielson's methods, an isovalue change needs to compute a new configuration.
The proposed method has been tested by using a cell with configurations for all the cases of improved marching cubes. In this section we show the most representative ones. Every example is shown using different isovalues on the basis of the same cell, where you can see how the isosurface changes continuously when the isovalue is slightly changed.
Example 1 (figure 10): A configuration with 1 component on the basis of the simpliest tessellation.
Example 2 (figure 13): By modifying the isovalue from 26 to 23 we can see how an isosurface of 2 components changes to an isosurface of 1 component (with a tunnel between them) around the body saddle point which is inside this ambiguous cell.
Example 3 (figure 14): By changing the isovalue from 40 to 26 we can see how the isosurface changes continuously from 2 components to 1 component through the face saddle point on the bottom face.
Example 4 (figure 15): By modifying the isovalue from 36 to 50 the isosurface of 2 components on two opposite ambiguous faces changes continuously to an isosurface of 1 component and again, to an isosurface of 2 components.
Example 5 (figure 16): This example shows a cell with all its faces ambiguous but without body saddle points. The isovalue changes from 43 to 63.
Example 6 (figure 18): This example shows a cell with 2 body saddle points, by changing the isovalue from 42 to 58 we can see how the isosurface forms a tunnel around each body saddle point depending on the particular isovalue.
The vertex values, in addition to the number of tetrahedra, of every example, are shown in table 1.
Table 1. Example's vertex values and number of tetrahedra (last row).
V. 
Ex.1 
Ex.2 
Ex.3 
Ex.4 
Ex.5 
Ex.6 
V_{0} 
0 
35 
50 
100 
78 
33 
V_{1} 
10 
12 
14 
71 
33 
100 
V_{2} 
20 
0 
14 
14 
33 
100 
V_{3} 
30 
38 
14 
14 
100 
0 
V_{4} 
40 
100 
0 
14 
0 
100 
V_{5} 
50 
15 
14 
0 
78 
0 
V_{6} 
60 
15 
100 
71 
78 
0 
V_{7} 
100 
15 
14 
71 
33 
67 

6 
14 
14 
16 
24 
30 
In general the method produces all Lopes' topologies correctly; in fact, all the topologies in figure 7 have been achieved using our method. As the EMC method by Lopes produces topollogicaly correct isosurfaces, our method also produces topollogicaly correct isosurfaces.
The method has also been tested with real data. Figure 19a shows a model of some teeth rendered using our method. It is possible to see that there are no holes on he surface. Figures 19b and 19c show a zoom of the image generated with our method and marching cubes respectively.
Figure 19. Visualization of a model of some teeth. (a) Visualization without zoom. (b) Zoomed visualization with the correct topology. (c) Zoomed visualization without the correct topology; some holes appear.
Our method produces similar results to Lopes [ LB03 ]. However, every time the isovalue changes, Lopes' method has to compute the new cell configuration (from a total of 31 distinct configurations), the new possible bishoulder points, and the new possible tangent points, and then the triangle configuration is formed and rendered.
Using Lopes method, some infomation can be precomputed just once, but other information depends on the isovalue, so this kind of information needs to be computed every time the isovalue changes. As can be seen in his work, the computation of bishoulder points needs to compute a minimum of 2 face shoulder points (or more if a more accurate approximation is needed) so several square roots need to be computed. The computation of tangent points also needs to compute square roots because three quadratic equations have to be solved even though they have the same discriminant. So, every time the isovalue changes, several costly computations have to be performed. That is, it is a time expensive method for interaction.
Our method also needs to compute square roots (just one per cell) to determine which tessellation must be carried out, however, the particular tessellation is computed just once for every cell because our tessellation is isovalue independent. As the tessellation does not depend on the classification of the vertices, once the tessellation has been carried out, the volume is represented by a set of tetrahedra and it is visualized by the well known and fast marching tetrahedra. No new tessellation needs to be computed, unless some value is changed. In this case, only the cells that share that point will be tessellated again. An isovalue change does not require any new tessellation.
Both methods are topologically valid but our method needs less computations when the isovalue changes as it is isovalue independent. So, it is faster with, allowing interactive visualization.
Moreover, our method also produces continuous changes in the isosurface when the isovalue changes continuously. It is so, because the saddle points, which are points where one topology changes into another one, are vertices of the tetrahedra. So, when the isovalue is changed, the isosurface is moved towards or from saddle points where its topology is changed
Hence, our three goals of section 4 are achieved:

The tessellation only depends on the cell vertex' values. It is independent of the isovalue.

The isosurfaces that are obtained are topollogically equivalents to the isosurfaces that are obtained by Lopes' method, so they are topollogically corrects.

The movement of the isosurface is continuous when the isovalue is moved continuously.
The drawback of our method is the number of tetrahedra that is produced, however acording to several test, the most frequent case is the case 00, which appears more than a 52%. The case 11 (which includes the case 10) appears a 27%, the case 21 (plus the case 20) appears a 13% and so on until the worst cases (cases n2) which, by taking into account all together, appear less than 1%. The tests have been carried out building cells within random values. Table 2 shows the results.
Case 
Frequency 
Case 00 
52.77% 
Cases 10 and 11 
27.22% 
Cases 20 and 21 
13.23% 
Cases 30 and 31 
3,39% 
Cases 40 and 41 
1.98% 
Cases 50 and 51 
0.00% 
Cases 60 and 61 
0.98% 
Cases 02, 12, 22 and 32 
0.00% 
Case 42 
0 000011% 
Case 52 
0.00% 
Case 62 
0.43% 
As future work, we want to analyse the possibility of grouping the tetrahedra in order to reduce their global amount and to allow a multiresolution representation of the volume.
This work has been funded by the Spanish Government and by ERDF funds under project TIN200406326C0302. It has also been unded by Andalusian Government under excellence project TIC401.
[CGMS00] Reconstruction of topologically correct and adaptive trilineal surfaces, Computer and Graphics, (2000), no. 3, 399—418, issn 00978493.
[Che95] Technical Report CN/9517, CERN, 1995 , http://wwwinfo.cern.ch/asdoc/psdir/mc.ps.gz, Last visited March 20th, 2008. Marching cubes 33: Construccion of Topologically Correct Isosurfaces,
[CL03] Progressive simplification of tetrahedral meshes preserving all isosurface topologies, Computer Graphics Forum, (2003), no. 3, 493—504, issn 01677055.
[CP98] A new tetrahedral tesselation scheme for isosurface generation, Computers & Graphics, (1998), no. 1, 83—90, issn 00978493.
[GH95] Exploiting triangulated surface extraction using tetrahedral decomposition, IEEE Transactions on Visualization and Computer Graphics, (1995), no. 4, 328—342, issn 10772626.
[GP00] Topology preserving and controlled topology simplifying multiresolution isosurface extraction, Proceedings of the IEEE Conference on Visualization'00, 2000, 259—266, Salt Lake City, Utah, USA, isbn 158113309X.
[LB03] Improving the robustness and accuracy of the marching cubes algorithm for isosurfacing, IEEE Transaction on Visualization and Computer Graphics, (2003), no. 1, 16—29, issn 10772626.
[LC87] Marching cubes: A high resolution 3d surface construction algorithm, ACM Computer Graphics, (1987), no. 4, 163—169, issn 00978930.
[Lev90] Efficient ray tracing of volume data, ACM Transactions on Graphics, (1990), no. 3, 245—261, issn 07300301.
[Lop99] Accuracy in Scientific Visualization, Phd Thesis, University of Leeds, United Kingdom, 1999.
[MSS94] A modified lookup table for implicit disambiguation of marching cubes, The Visual Computer, (1994), no. 6, 353—355, issn 01782789.
[NH91] The asymptotic decider: resolving the ambiguity in marching cubes, Proceedings of the IEEE Conference on Visualization'91, pp. 83—91, 1991, San Diego, California, USA, isbn 0818622458.
[Nie03] On marching cubes, IEEE Transactions on visualization and computer graphics, (2003), no. 3, 283—297, issn 10772626.
[NY06] A survey of the marching cubes algorithm, Computers & Graphics, (2006), no. 5, 854—879, issn 00978493.
[PT90] Surface mapping brain function on 3d models, IEEE Computer Graphics & Applications, (1990), no. 5, 33—41, issn 02721716.
[TPG98] Regularised marching tetrahedra: Improved isosurface extraction, Computers & Graphics, (1999), no. 4, 583—598, issn 00978493.
[WG90] Topological considerations in isosurface generation. extended abstract, ACM Computer Graphics, (1990), no. 5, 79—86, issn 00978930.
[ZCT95] An elaborate ambiguity detection method for constructing isosurfaces within tetrahedral meshes, Computers & Graphics, (1995), no. 3, 355—364, issn 00978493.
Fulltext ¶
 Volltext als PDF ( Size 3.4 MB )
License ¶
Any party may pass on this Work by electronic means and make it available for download under the terms and conditions of the Digital Peer Publishing License. The text of the license may be accessed and retrieved at http://www.dipp.nrw.de/lizenzen/dppl/dppl/DPPL_v2_en_062004.html.
Recommended citation ¶
Francisco Velasco, Juan Carlos Torres, Alejandro Léon, and Francisco Soler, Adaptive Cube Tessellation for Topologically Correct Isosurfaces. JVRB  Journal of Virtual Reality and Broadcasting, 5(2008), no. 3. (urn:nbn:de:0009613098)
Please provide the exact URL and date of your last visit when citing this article.