25 #ifndef WMARCHINGCUBESALGORITHM_H
26 #define WMARCHINGCUBESALGORITHM_H
30 #include "../../common/math/WMatrix.h"
31 #include "../../common/WProgressCombiner.h"
32 #include "../WTriangleMesh.h"
33 #include "marchingCubesCaseTables.h"
47 typedef std::map< unsigned int, WPointXYZId > ID2WPointXYZId;
57 typedef std::vector<WMCTriangle> WMCTriangleVECTOR;
113 template<
typename T >
114 boost::shared_ptr< WTriangleMesh >
generateSurface(
size_t nbCoordsX,
size_t nbCoordsY,
size_t nbCoordsZ,
116 const std::vector< T >* vals,
118 boost::shared_ptr< WProgressCombiner > mainProgress );
135 unsigned int nX,
unsigned int nY,
unsigned int nZ,
unsigned int nEdgeNo );
152 WPointXYZId interpolate(
double fX1,
double fY1,
double fZ1,
double fX2,
double fY2,
double fZ2,
double tVal1,
double tVal2 );
164 int getEdgeID(
unsigned int nX,
unsigned int nY,
unsigned int nZ,
unsigned int nEdgeNo );
174 unsigned int getVertexID(
unsigned int nX,
unsigned int nY,
unsigned int nZ );
191 const std::vector< T >* vals,
193 boost::shared_ptr< WProgressCombiner > mainProgress )
195 WAssert( vals,
"No value set provided." );
212 unsigned int nPointsInSlice = nX * nY;
214 boost::shared_ptr< WProgress > progress = boost::shared_ptr< WProgress >(
new WProgress(
"Marching Cubes",
m_nCellsZ ) );
215 mainProgress->addSubProgress( progress );
217 for(
unsigned int z = 0; z <
m_nCellsZ; z++ )
220 for(
unsigned int y = 0; y <
m_nCellsY; y++ )
222 for(
unsigned int x = 0; x <
m_nCellsX; x++ )
226 unsigned int tableIndex = 0;
227 if( ( *vals )[ z * nPointsInSlice + y * nX + x ] <
m_tIsoLevel )
229 if( ( *vals )[ z * nPointsInSlice + ( y + 1 ) * nX + x ] <
m_tIsoLevel )
231 if( ( *vals )[ z * nPointsInSlice + ( y + 1 ) * nX + ( x + 1 ) ] <
m_tIsoLevel )
233 if( ( *vals )[ z * nPointsInSlice + y * nX + ( x + 1 ) ] <
m_tIsoLevel )
235 if( ( *vals )[ ( z + 1 ) * nPointsInSlice + y * nX + x ] <
m_tIsoLevel )
237 if( ( *vals )[ ( z + 1 ) * nPointsInSlice + ( y + 1 ) * nX + x ] <
m_tIsoLevel )
239 if( ( *vals )[ ( z + 1 ) * nPointsInSlice + ( y + 1 ) * nX + ( x + 1 ) ] <
m_tIsoLevel )
241 if( ( *vals )[ ( z + 1 ) * nPointsInSlice + y * nX + ( x + 1 ) ] <
m_tIsoLevel )
246 if( wMarchingCubesCaseTables::edgeTable[tableIndex] != 0 )
248 if( wMarchingCubesCaseTables::edgeTable[tableIndex] & 8 )
251 unsigned int id =
getEdgeID( x, y, z, 3 );
254 if( wMarchingCubesCaseTables::edgeTable[tableIndex] & 1 )
257 unsigned int id =
getEdgeID( x, y, z, 0 );
260 if( wMarchingCubesCaseTables::edgeTable[tableIndex] & 256 )
263 unsigned int id =
getEdgeID( x, y, z, 8 );
267 if( x == m_nCellsX - 1 )
269 if( wMarchingCubesCaseTables::edgeTable[tableIndex] & 4 )
272 unsigned int id =
getEdgeID( x, y, z, 2 );
275 if( wMarchingCubesCaseTables::edgeTable[tableIndex] & 2048 )
278 unsigned int id =
getEdgeID( x, y, z, 11 );
282 if( y == m_nCellsY - 1 )
284 if( wMarchingCubesCaseTables::edgeTable[tableIndex] & 2 )
287 unsigned int id =
getEdgeID( x, y, z, 1 );
290 if( wMarchingCubesCaseTables::edgeTable[tableIndex] & 512 )
293 unsigned int id =
getEdgeID( x, y, z, 9 );
297 if( z == m_nCellsZ - 1 )
299 if( wMarchingCubesCaseTables::edgeTable[tableIndex] & 16 )
302 unsigned int id =
getEdgeID( x, y, z, 4 );
305 if( wMarchingCubesCaseTables::edgeTable[tableIndex] & 128 )
308 unsigned int id =
getEdgeID( x, y, z, 7 );
312 if( ( x == m_nCellsX - 1 ) && ( y == m_nCellsY - 1 ) )
313 if( wMarchingCubesCaseTables::edgeTable[tableIndex] & 1024 )
316 unsigned int id =
getEdgeID( x, y, z, 10 );
319 if( ( x == m_nCellsX - 1 ) && ( z == m_nCellsZ - 1 ) )
320 if( wMarchingCubesCaseTables::edgeTable[tableIndex] & 64 )
323 unsigned int id =
getEdgeID( x, y, z, 6 );
326 if( ( y == m_nCellsY - 1 ) && ( z == m_nCellsZ - 1 ) )
327 if( wMarchingCubesCaseTables::edgeTable[tableIndex] & 32 )
330 unsigned int id =
getEdgeID( x, y, z, 5 );
334 for(
int i = 0; wMarchingCubesCaseTables::triTable[tableIndex][i] != -1; i += 3 )
337 unsigned int pointID0, pointID1, pointID2;
338 pointID0 =
getEdgeID( x, y, z, wMarchingCubesCaseTables::triTable[tableIndex][i] );
339 pointID1 =
getEdgeID( x, y, z, wMarchingCubesCaseTables::triTable[tableIndex][i + 1] );
340 pointID2 =
getEdgeID( x, y, z, wMarchingCubesCaseTables::triTable[tableIndex][i + 2] );
341 triangle.
pointID[0] = pointID0;
342 triangle.
pointID[1] = pointID1;
343 triangle.
pointID[2] = pointID2;
351 unsigned int nextID = 0;
359 mapIterator->second.y / nbCoordsY,
360 mapIterator->second.z / nbCoordsZ );
363 WPosition pos =
WPosition( mapIterator->second.x, mapIterator->second.y, mapIterator->second.z );
365 std::vector< double > resultPos4D( 4 );
371 ( *mapIterator ).second.newID = nextID;
372 triMesh->addVertex( resultPos4D[0] / resultPos4D[3], resultPos4D[1] / resultPos4D[3], resultPos4D[2] / resultPos4D[3] );
373 triMesh->addTextureCoordinate( texCoord );
382 for(
unsigned int i = 0; i < 3; i++ )
384 unsigned int newID =
m_idToVertices[( *vecIterator ).pointID[i]].newID;
385 ( *vecIterator ).pointID[i] = newID;
387 triMesh->addTriangle( ( *vecIterator ).pointID[0], ( *vecIterator ).pointID[1], ( *vecIterator ).pointID[2] );
396 unsigned int nX,
unsigned int nY,
unsigned int nZ,
397 unsigned int nEdgeNo )
407 unsigned int v1x = nX;
408 unsigned int v1y = nY;
409 unsigned int v1z = nZ;
411 unsigned int v2x = nX;
412 unsigned int v2y = nY;
413 unsigned int v2z = nZ;
487 double val1 = ( *vals )[ v1z * nPointsInSlice + v1y * (
m_nCellsX + 1 ) + v1x ];
488 double val2 = ( *vals )[ v2z * nPointsInSlice + v2y * (
m_nCellsX + 1 ) + v2x ];
491 intersection.newID = 0;
495 #endif // WMARCHINGCUBESALGORITHM_H