iRoCS Toolbox  1.1.0
MarchingCubes.hh
Go to the documentation of this file.
1 /**************************************************************************
2  *
3  * Copyright (C) 2015 Thorsten Falk
4  *
5  * Image Analysis Lab, University of Freiburg, Germany
6  *
7  * This program is free software; you can redistribute it and/or modify
8  * it under the terms of the GNU General Public License as published by
9  * the Free Software Foundation; either version 3 of the License, or
10  * (at your option) any later version.
11  *
12  * This program is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  * GNU General Public License for more details.
16  *
17  * You should have received a copy of the GNU General Public License
18  * along with this program; if not, write to the Free Software Foundation,
19  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
20  *
21  **************************************************************************/
22 
23 /*======================================================================*/
29 /*======================================================================*/
30 
31 #ifndef ATBMARCHINGCUBES_HH
32 #define ATBMARCHINGCUBES_HH
33 
34 #ifdef HAVE_CONFIG_H
35 #include <config.hh>
36 #endif
37 
38 #include <blitz/array.h>
39 
40 #include <vector>
41 #include <set>
42 
43 #include "Array.hh"
44 #include "SurfaceGeometry.hh"
45 #include "Neighborhood.hh"
46 
47 namespace atb
48 {
49 
50 /*======================================================================*/
62 /*======================================================================*/
64  {
65 
66  public:
67 
68 /*======================================================================*/
72 /*======================================================================*/
73  typedef blitz::TinyVector<SurfaceGeometry::VertexT,3> Triangle;
74 
75  private:
76 
77  struct GridCell
78  {
79  blitz::TinyVector<SurfaceGeometry::VertexT,8> cornerCoordinates;
80  blitz::TinyVector<double,8> cornerValues;
81  };
82 
83  static int edgeTable[256];
84  static int triangleTable[256][16];
85 
86  static int computeTrianglesForGridCell(
87  GridCell const &gridCell, std::vector<Triangle> &triangles,
88  double isoLevel = 0.0);
89 
90  static void computeTrianglesForGridCell(
91  GridCell const &gridCell, blitz::TinyVector<int,12> &cellIndices,
92  SurfaceGeometry &surface, double isoLevel = 0.0);
93 
94  static SurfaceGeometry::VertexT getInterpolatedVertex(
95  double isoLevel, SurfaceGeometry::VertexT const &p1,
96  SurfaceGeometry::VertexT const &p2, double val1, double val2);
97 
98  public:
99 
100 /*======================================================================*/
119 /*======================================================================*/
120  template<typename DataT>
121  static void triangulate(
122  Array<DataT,3> const &data, SurfaceGeometry &surface,
123  double isoLevel = 0.0, double simplifyTolerance = 0.0);
124 
125 /*======================================================================*/
145 /*======================================================================*/
146  template<typename DataT>
147  static void triangulate(
148  blitz::Array<DataT,3> const &data,
149  blitz::TinyVector<double,3> const &elementSizeUm,
150  SurfaceGeometry &surface, double isoLevel = 0.0,
151  double simplifyTolerance = 0.0);
152 
153 /*======================================================================*/
164 /*======================================================================*/
165  static void simplifyMesh(
166  SurfaceGeometry &surface, double simplifyTolerance);
167 
168  private:
169 
170  static double _computeEdgeRemovalCost(
171  std::pair<size_t,size_t> const &edge, SurfaceGeometry const &surface,
172  std::vector< std::set<size_t> > const &triangles);
173 
174  static void _checkAndFixMesh(SurfaceGeometry &surface);
175 
176  };
177 
178 }
179 
180 #include "MarchingCubes.icc"
181 
182 #endif
Array class derived from blitz++ Arrays for handling microscopic datasets with associated element siz...
The Array class is an extension to the blitz++ Array class providing additional parameters element si...
Definition: Array.hh:85
blitz::TinyVector< float, 3 > VertexT
VertexT is the type used to store a 3D vertex position.
Storage and rendering of Triangulated surfaces.
Neighborhoods for local operators.
static void triangulate(Array< DataT, 3 > const &data, SurfaceGeometry &surface, double isoLevel=0.0, double simplifyTolerance=0.0)
The triangulate function computes a vertex mesh from a given levelset Array.
blitz::TinyVector< SurfaceGeometry::VertexT, 3 > Triangle
A triangle in 3-D space.
static void simplifyMesh(SurfaceGeometry &surface, double simplifyTolerance)
Simplify the given mesh by merging vertices as long as the given tolerance is not exceeded...
The SurfaceGeometry struct provides data structures required for the storage of triangulated surface ...
The MarchingCubes class provides the Marching Cubes Algorithm for isosurface tesselation of real scal...