OpenWalnut 1.2.5
WMarchingLegoAlgorithm.cpp
00001 //---------------------------------------------------------------------------
00002 //
00003 // Project: OpenWalnut ( http://www.openwalnut.org )
00004 //
00005 // Copyright 2009 OpenWalnut Community, BSV@Uni-Leipzig and CNCF@MPI-CBS
00006 // For more information see http://www.openwalnut.org/copying
00007 //
00008 // This file is part of OpenWalnut.
00009 //
00010 // OpenWalnut is free software: you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as published by
00012 // the Free Software Foundation, either version 3 of the License, or
00013 // (at your option) any later version.
00014 //
00015 // OpenWalnut is distributed in the hope that it will be useful,
00016 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00018 // GNU Lesser General Public License for more details.
00019 //
00020 // You should have received a copy of the GNU Lesser General Public License
00021 // along with OpenWalnut. If not, see <http://www.gnu.org/licenses/>.
00022 //
00023 //---------------------------------------------------------------------------
00024 
00025 #include "WMarchingLegoAlgorithm.h"
00026 
00027 WMarchingLegoAlgorithm::WMarchingLegoAlgorithm()
00028     : m_matrix( 4, 4 )
00029 {
00030 }
00031 
00032 WMarchingLegoAlgorithm::~WMarchingLegoAlgorithm()
00033 {
00034 }
00035 
00036 void WMarchingLegoAlgorithm::addSurface( size_t x, size_t y, size_t z, size_t surface )
00037 {
00038     WMLPointXYZId pt1;
00039     WMLPointXYZId pt2;
00040     WMLPointXYZId pt3;
00041     WMLPointXYZId pt4;
00042 
00043     pt1.newID = 0;
00044     pt2.newID = 0;
00045     pt3.newID = 0;
00046     pt4.newID = 0;
00047 
00048     switch ( surface )
00049     {
00050         case 1:
00051         {
00052             pt1.x = x;
00053             pt1.y = y;
00054             pt1.z = z;
00055             unsigned int id1 = getVertexID( pt1.x , pt1.y, pt1.z );
00056             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id1, pt1 ) );
00057 
00058             pt2.x = x;
00059             pt2.y = y + 1;
00060             pt2.z = z;
00061             unsigned int id2 = getVertexID( pt2.x , pt2.y, pt2.z );
00062             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id2, pt2 ) );
00063 
00064             pt3.x = x;
00065             pt3.y = y + 1;
00066             pt3.z = z + 1;
00067             unsigned int id3 = getVertexID( pt3.x , pt3.y, pt3.z );
00068             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id3, pt3 ) );
00069 
00070             pt4.x = x;
00071             pt4.y = y;
00072             pt4.z = z + 1;
00073             unsigned int id4 = getVertexID( pt4.x , pt4.y, pt4.z );
00074             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id4, pt4 ) );
00075 
00076             WMLTriangle triangle1;
00077             triangle1.pointID[0] = id1;
00078             triangle1.pointID[1] = id2;
00079             triangle1.pointID[2] = id3;
00080             WMLTriangle triangle2;
00081             triangle2.pointID[0] = id3;
00082             triangle2.pointID[1] = id4;
00083             triangle2.pointID[2] = id1;
00084             m_trivecTriangles.push_back( triangle1 );
00085             m_trivecTriangles.push_back( triangle2 );
00086             break;
00087         }
00088         case 2:
00089         {
00090             pt1.x = x + 1;
00091             pt1.y = y;
00092             pt1.z = z;
00093             unsigned int id1 = getVertexID( pt1.x , pt1.y, pt1.z );
00094             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id1, pt1 ) );
00095 
00096             pt2.x = x + 1;
00097             pt2.y = y;
00098             pt2.z = z + 1;
00099             unsigned int id2 = getVertexID( pt2.x , pt2.y, pt2.z );
00100             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id2, pt2 ) );
00101 
00102             pt3.x = x + 1;
00103             pt3.y = y + 1;
00104             pt3.z = z + 1;
00105             unsigned int id3 = getVertexID( pt3.x , pt3.y, pt3.z );
00106             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id3, pt3 ) );
00107 
00108             pt4.x = x + 1;
00109             pt4.y = y + 1;
00110             pt4.z = z;
00111             unsigned int id4 = getVertexID( pt4.x , pt4.y, pt4.z );
00112             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id4, pt4 ) );
00113 
00114             WMLTriangle triangle1;
00115             triangle1.pointID[0] = id1;
00116             triangle1.pointID[1] = id2;
00117             triangle1.pointID[2] = id3;
00118             WMLTriangle triangle2;
00119             triangle2.pointID[0] = id3;
00120             triangle2.pointID[1] = id4;
00121             triangle2.pointID[2] = id1;
00122             m_trivecTriangles.push_back( triangle1 );
00123             m_trivecTriangles.push_back( triangle2 );
00124             break;
00125         }
00126         case 3:
00127         {
00128             pt1.x = x;
00129             pt1.y = y;
00130             pt1.z = z;
00131             unsigned int id1 = getVertexID( pt1.x , pt1.y, pt1.z );
00132             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id1, pt1 ) );
00133 
00134             pt2.x = x;
00135             pt2.y = y;
00136             pt2.z = z + 1;
00137             unsigned int id2 = getVertexID( pt2.x , pt2.y, pt2.z );
00138             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id2, pt2 ) );
00139 
00140             pt3.x = x + 1;
00141             pt3.y = y;
00142             pt3.z = z + 1;
00143             unsigned int id3 = getVertexID( pt3.x , pt3.y, pt3.z );
00144             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id3, pt3 ) );
00145 
00146             pt4.x = x + 1;
00147             pt4.y = y;
00148             pt4.z = z;
00149             unsigned int id4 = getVertexID( pt4.x , pt4.y, pt4.z );
00150             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id4, pt4 ) );
00151 
00152             WMLTriangle triangle1;
00153             triangle1.pointID[0] = id1;
00154             triangle1.pointID[1] = id2;
00155             triangle1.pointID[2] = id3;
00156             WMLTriangle triangle2;
00157             triangle2.pointID[0] = id3;
00158             triangle2.pointID[1] = id4;
00159             triangle2.pointID[2] = id1;
00160             m_trivecTriangles.push_back( triangle1 );
00161             m_trivecTriangles.push_back( triangle2 );
00162             break;
00163         }
00164         case 4:
00165         {
00166             pt1.x = x;
00167             pt1.y = y + 1;
00168             pt1.z = z;
00169             unsigned int id1 = getVertexID( pt1.x , pt1.y, pt1.z );
00170             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id1, pt1 ) );
00171 
00172             pt2.x = x + 1;
00173             pt2.y = y + 1;
00174             pt2.z = z;
00175             unsigned int id2 = getVertexID( pt2.x , pt2.y, pt2.z );
00176             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id2, pt2 ) );
00177 
00178             pt3.x = x + 1;
00179             pt3.y = y + 1;
00180             pt3.z = z + 1;
00181             unsigned int id3 = getVertexID( pt3.x , pt3.y, pt3.z );
00182             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id3, pt3 ) );
00183 
00184             pt4.x = x;
00185             pt4.y = y + 1;
00186             pt4.z = z + 1;
00187             unsigned int id4 = getVertexID( pt4.x , pt4.y, pt4.z );
00188             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id4, pt4 ) );
00189 
00190             WMLTriangle triangle1;
00191             triangle1.pointID[0] = id1;
00192             triangle1.pointID[1] = id2;
00193             triangle1.pointID[2] = id3;
00194             WMLTriangle triangle2;
00195             triangle2.pointID[0] = id3;
00196             triangle2.pointID[1] = id4;
00197             triangle2.pointID[2] = id1;
00198             m_trivecTriangles.push_back( triangle1 );
00199             m_trivecTriangles.push_back( triangle2 );
00200             break;
00201         }
00202         case 5:
00203         {
00204             pt1.x = x;
00205             pt1.y = y;
00206             pt1.z = z;
00207             unsigned int id1 = getVertexID( pt1.x , pt1.y, pt1.z );
00208             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id1, pt1 ) );
00209 
00210             pt2.x = x + 1;
00211             pt2.y = y;
00212             pt2.z = z;
00213             unsigned int id2 = getVertexID( pt2.x , pt2.y, pt2.z );
00214             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id2, pt2 ) );
00215 
00216             pt3.x = x + 1;
00217             pt3.y = y + 1;
00218             pt3.z = z;
00219             unsigned int id3 = getVertexID( pt3.x , pt3.y, pt3.z );
00220             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id3, pt3 ) );
00221 
00222             pt4.x = x;
00223             pt4.y = y + 1;
00224             pt4.z = z;
00225             unsigned int id4 = getVertexID( pt4.x , pt4.y, pt4.z );
00226             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id4, pt4 ) );
00227 
00228             WMLTriangle triangle1;
00229             triangle1.pointID[0] = id1;
00230             triangle1.pointID[1] = id2;
00231             triangle1.pointID[2] = id3;
00232             WMLTriangle triangle2;
00233             triangle2.pointID[0] = id3;
00234             triangle2.pointID[1] = id4;
00235             triangle2.pointID[2] = id1;
00236             m_trivecTriangles.push_back( triangle1 );
00237             m_trivecTriangles.push_back( triangle2 );
00238             break;
00239         }
00240         case 6:
00241         {
00242             pt1.x = x;
00243             pt1.y = y;
00244             pt1.z = z + 1;
00245             unsigned int id1 = getVertexID( pt1.x , pt1.y, pt1.z );
00246             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id1, pt1 ) );
00247 
00248             pt2.x = x;
00249             pt2.y = y + 1;
00250             pt2.z = z + 1;
00251             unsigned int id2 = getVertexID( pt2.x , pt2.y, pt2.z );
00252             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id2, pt2 ) );
00253 
00254             pt3.x = x + 1;
00255             pt3.y = y + 1;
00256             pt3.z = z + 1;
00257             unsigned int id3 = getVertexID( pt3.x , pt3.y, pt3.z );
00258             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id3, pt3 ) );
00259 
00260             pt4.x = x + 1;
00261             pt4.y = y;
00262             pt4.z = z + 1;
00263             unsigned int id4 = getVertexID( pt4.x , pt4.y, pt4.z );
00264             m_idToVertices.insert( ID2WMLPointXYZId::value_type( id4, pt4 ) );
00265 
00266             WMLTriangle triangle1;
00267             triangle1.pointID[0] = id1;
00268             triangle1.pointID[1] = id2;
00269             triangle1.pointID[2] = id3;
00270             WMLTriangle triangle2;
00271             triangle2.pointID[0] = id3;
00272             triangle2.pointID[1] = id4;
00273             triangle2.pointID[2] = id1;
00274             m_trivecTriangles.push_back( triangle1 );
00275             m_trivecTriangles.push_back( triangle2 );
00276             break;
00277         }
00278         default:
00279             break;
00280     }
00281 }
00282 
00283 size_t WMarchingLegoAlgorithm::getVertexID( size_t nX, size_t nY, size_t nZ )
00284 {
00285     return nZ * ( m_nCellsY + 1 ) * ( m_nCellsX + 1) + nY * ( m_nCellsX + 1 ) + nX;
00286 }
00287 
00288 boost::shared_ptr<WTriangleMesh> WMarchingLegoAlgorithm::genSurfaceOneValue( size_t nbCoordsX, size_t nbCoordsY, size_t nbCoordsZ,
00289                                                                                                  const WMatrix< double >& mat,
00290                                                                                                  const std::vector< size_t >* vals,
00291                                                                                                  size_t isoValue )
00292 {
00293     WAssert( vals, "No value set provided." );
00294 
00295     m_idToVertices.clear();
00296     m_trivecTriangles.clear();
00297 
00298     m_nCellsX = nbCoordsX - 1;
00299     m_nCellsY = nbCoordsY - 1;
00300     m_nCellsZ = nbCoordsZ - 1;
00301 
00302     m_matrix = mat;
00303 
00304     size_t nX = nbCoordsX;
00305     size_t nY = nbCoordsY;
00306 
00307     size_t nPointsInSlice = nX * nY;
00308 
00309     // Generate isosurface.
00310     for( size_t z = 0; z < m_nCellsZ; z++ )
00311     {
00312         for( size_t y = 0; y < m_nCellsY; y++ )
00313         {
00314             for( size_t x = 0; x < m_nCellsX; x++ )
00315             {
00316                 if( ( *vals )[ z * nPointsInSlice + y * nX + x ] != isoValue )
00317                 {
00318                     continue;
00319                 }
00320 
00321                 if( x > 0 && ( ( *vals )[ z * nPointsInSlice + y * nX + x - 1 ] != isoValue ) )
00322                 {
00323                     addSurface( x, y, z, 1 );
00324                 }
00325                 if( x < m_nCellsX - 1 && ( ( *vals )[ z * nPointsInSlice + y * nX + x + 1 ] != isoValue ) )
00326                 {
00327                     addSurface( x, y, z, 2 );
00328                 }
00329 
00330                 if( y > 0 && ( ( *vals )[ z * nPointsInSlice + ( y - 1 ) * nX + x ] != isoValue ) )
00331                 {
00332                     addSurface( x, y, z, 3 );
00333                 }
00334 
00335                 if( y < m_nCellsY - 1 && ( ( *vals )[ z * nPointsInSlice + ( y + 1 ) * nX + x ] != isoValue ) )
00336                 {
00337                     addSurface( x, y, z, 4 );
00338                 }
00339 
00340                 if( z > 0 && ( ( *vals )[ ( z - 1 ) * nPointsInSlice + y * nX + x ] != isoValue ) )
00341                 {
00342                     addSurface( x, y, z, 5 );
00343                 }
00344 
00345                 if( z < m_nCellsZ - 1 && ( ( *vals )[ ( z + 1 ) * nPointsInSlice + y * nX + x ] != isoValue ) )
00346                 {
00347                     addSurface( x, y, z, 6 );
00348                 }
00349 
00350                 if( x == 0 )
00351                 {
00352                     addSurface( x, y, z, 1 );
00353                 }
00354                 if( x == m_nCellsX - 1 )
00355                 {
00356                     addSurface( x, y, z, 2 );
00357                 }
00358 
00359                 if( y == 0 )
00360                 {
00361                     addSurface( x, y, z, 3 );
00362                 }
00363 
00364                 if( y == m_nCellsY - 1 )
00365                 {
00366                     addSurface( x, y, z, 4 );
00367                 }
00368 
00369                 if( z == 0 )
00370                 {
00371                     addSurface( x, y, z, 5 );
00372                 }
00373 
00374                 if( z == m_nCellsZ - 1 )
00375                 {
00376                     addSurface( x, y, z, 6 );
00377                 }
00378             }
00379         }
00380     }
00381     unsigned int nextID = 0;
00382     boost::shared_ptr< WTriangleMesh > triMesh( new WTriangleMesh( m_idToVertices.size(), m_trivecTriangles.size() ) );
00383 
00384     // Rename vertices.
00385     ID2WMLPointXYZId::iterator mapIterator = m_idToVertices.begin();
00386     while( mapIterator != m_idToVertices.end() )
00387     {
00388         WPosition texCoord = WPosition( mapIterator->second.x / nbCoordsX,
00389                                                       mapIterator->second.y / nbCoordsY,
00390                                                       mapIterator->second.z / nbCoordsZ );
00391 
00392         // transform from grid coordinate system to world coordinates
00393         WPosition pos = WPosition( mapIterator->second.x, mapIterator->second.y, mapIterator->second.z );
00394 
00395         std::vector< double > resultPos4D( 4 );
00396         resultPos4D[0] = m_matrix( 0, 0 ) * pos[0] + m_matrix( 0, 1 ) * pos[1] + m_matrix( 0, 2 ) * pos[2] + m_matrix( 0, 3 ) * 1;
00397         resultPos4D[1] = m_matrix( 1, 0 ) * pos[0] + m_matrix( 1, 1 ) * pos[1] + m_matrix( 1, 2 ) * pos[2] + m_matrix( 1, 3 ) * 1;
00398         resultPos4D[2] = m_matrix( 2, 0 ) * pos[0] + m_matrix( 2, 1 ) * pos[1] + m_matrix( 2, 2 ) * pos[2] + m_matrix( 2, 3 ) * 1;
00399         resultPos4D[3] = m_matrix( 3, 0 ) * pos[0] + m_matrix( 3, 1 ) * pos[1] + m_matrix( 3, 2 ) * pos[2] + m_matrix( 3, 3 ) * 1;
00400 
00401         ( *mapIterator ).second.newID = nextID;
00402         triMesh->addVertex( resultPos4D[0] / resultPos4D[3],
00403                             resultPos4D[1] / resultPos4D[3],
00404                             resultPos4D[2] / resultPos4D[3] );
00405         triMesh->addTextureCoordinate( texCoord );
00406         nextID++;
00407         mapIterator++;
00408     }
00409 
00410     // Now rename triangles.
00411     WMLTriangleVECTOR::iterator vecIterator = m_trivecTriangles.begin();
00412     while( vecIterator != m_trivecTriangles.end() )
00413     {
00414         for( unsigned int i = 0; i < 3; i++ )
00415         {
00416             unsigned int newID = m_idToVertices[( *vecIterator ).pointID[i]].newID;
00417             ( *vecIterator ).pointID[i] = newID;
00418         }
00419         triMesh->addTriangle( ( *vecIterator ).pointID[0], ( *vecIterator ).pointID[1], ( *vecIterator ).pointID[2] );
00420         vecIterator++;
00421     }
00422     return triMesh;
00423 }
 All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Friends