DGtal  0.9.4beta
topology/volToOFF.cpp

Marching-cube like surface extracted using the combinatorial manifold structure of digital surfaces.

See also
Application to export surface in OFF format

On the lobser.vol volume, volToOFF.cpp builds an OFF surface of 155068 vertices, 154910 faces, 310136 edges in 3879ms+1646ms.

* # Commands
* $ ./examples/topology/volToOff ../examples/samples/cat10.vol 1 255 0
* $ ./examples/topology/volToOff ../examples/samples/lobster.vol 50 255 0
* 
digital-surface-mc-cat10.png
Marching-cube surface of cat10.vol file.
digital-surface-mc-lobster.png
Marching-cube surface of lobster.vol file.
#include <iostream>
#include "DGtal/base/Common.h"
#include "DGtal/io/readers/VolReader.h"
#include "DGtal/helpers/StdDefs.h"
#include "DGtal/topology/helpers/Surfaces.h"
#include "DGtal/topology/DigitalSurface.h"
#include "DGtal/topology/SetOfSurfels.h"
#include "DGtal/images/ImageSelector.h"
#include "DGtal/images/imagesSetsUtils/SetFromImage.h"
using namespace std;
using namespace DGtal;
using namespace Z3i;
void usage( int, char** argv )
{
std::cerr << "Usage: " << argv[ 0 ] << " <fileName.vol> <minT> <maxT> <int=0|ext=1>" << std::endl;
std::cerr << "\t - displays the boundary of the shape stored in vol file <fileName.vol> as an OFF geomview surface file. It is a kind of marching-cube surface, defined by duality with respect to the digital surface." << std::endl;
std::cerr << "\t - voxel v belongs to the shape iff its value I(v) follows minT <= I(v) <= maxT." << std::endl;
std::cerr << "\t - 0: interior adjacency, 1: exterior adjacency." << std::endl;
}
int main( int argc, char** argv )
{
if ( argc < 5 )
{
usage( argc, argv );
return 1;
}
std::string inputFilename = argv[ 1 ];
unsigned int minThreshold = atoi( argv[ 2 ] );
unsigned int maxThreshold = atoi( argv[ 3 ] );
bool intAdjacency = atoi( argv[ 4 ] ) == 0;
trace.beginBlock( "Reading vol file into an image." );
Image image = VolReader<Image>::importVol(inputFilename);
DigitalSet set3d (image.domain());
SetFromImage<DigitalSet>::append<Image>(set3d, image,
minThreshold, maxThreshold);
// Construct the Khalimsky space from the image domain
KSpace K;
bool space_ok = K.init( image.domain().lowerBound(),
image.domain().upperBound(), true );
if (!space_ok)
{
trace.error() << "Error in the Khamisky space construction."<<std::endl;
return 2;
}
typedef SurfelAdjacency<KSpace::dimension> MySurfelAdjacency;
MySurfelAdjacency surfAdj( intAdjacency ); // interior in all directions.
trace.beginBlock( "Extracting boundary by scanning the space. " );
typedef KSpace::SurfelSet SurfelSet;
typedef SetOfSurfels< KSpace, SurfelSet > MySetOfSurfels;
MySetOfSurfels theSetOfSurfels( K, surfAdj );
Surfaces<KSpace>::sMakeBoundary( theSetOfSurfels.surfelSet(),
K, set3d,
image.domain().lowerBound(),
image.domain().upperBound() );
MyDigitalSurface digSurf( theSetOfSurfels );
trace.info() << "Digital surface has " << digSurf.size() << " surfels."
<< std::endl;
trace.beginBlock( "Making OFF surface <marching-cube.off>. " );
ofstream out( "marching-cube.off" );
if ( out.good() )
digSurf.exportSurfaceAs3DOFF( out );
out.close();
return 0;
}