ITK  4.8.0
Insight Segmentation and Registration Toolkit
SphinxExamples/src/Segmentation/LevelSetsv4Visualization/VisualizeEvolvingDense2DLevelSetAsElevationMap/Code.cxx
int main( int argc, char* argv[] )
{
if( argc != 3)
{
std::cerr << "Missing Arguments" << std::endl;
std::cerr << argv[0] << std::endl;
std::cerr << "1- Input Image" << std::endl;
std::cerr << "2- Number of Iterations" << std::endl;
return EXIT_FAILURE;
}
// Image Dimension
const unsigned int Dimension = 2;
typedef unsigned char InputPixelType;
// Read input image (to be processed).
ReaderType::Pointer reader = ReaderType::New();
reader->SetFileName( argv[1] );
reader->Update();
InputImageType::Pointer input = reader->GetOutput();
int numberOfIterations = atoi( argv[2] );
typedef float LevelSetPixelType;
typedef itk::Image< LevelSetPixelType, Dimension > LevelSetImageType;
typedef LevelSetType::OutputType LevelSetOutputType;
typedef LevelSetType::OutputRealType LevelSetRealType;
// Generate a binary mask that will be used as initialization for the level
// set evolution.
BinaryImageType::Pointer binary = BinaryImageType::New();
binary->SetRegions( input->GetLargestPossibleRegion() );
binary->CopyInformation( input );
binary->Allocate();
BinaryImageType::RegionType region;
BinaryImageType::IndexType index;
BinaryImageType::SizeType size;
index.Fill( 5 );
size.Fill( 120 );
region.SetIndex( index );
region.SetSize( size );
InputIteratorType iIt( binary, region );
iIt.GoToBegin();
while( !iIt.IsAtEnd() )
{
++iIt;
}
// convert a binary mask to a level-set function
typedef itk::BinaryImageToLevelSetImageAdaptor< BinaryImageType,
LevelSetType > BinaryImageToLevelSetType;
BinaryImageToLevelSetType::Pointer adaptor = BinaryImageToLevelSetType::New();
adaptor->SetInputImage( binary );
adaptor->Initialize();
LevelSetType::Pointer levelSet = adaptor->GetLevelSet();
// The Heaviside function
HeavisideFunctionType::Pointer heaviside = HeavisideFunctionType::New();
heaviside->SetEpsilon( 1.5 );
// Create the level set container
LevelSetContainerType::Pointer levelSetContainer = LevelSetContainerType::New();
levelSetContainer->SetHeaviside( heaviside );
levelSetContainer->AddLevelSet( 0, levelSet );
// Create the terms.
//
// // Chan and Vese internal term
ChanAndVeseInternalTermType::Pointer cvInternalTerm = ChanAndVeseInternalTermType::New();
cvInternalTerm->SetInput( input );
cvInternalTerm->SetCoefficient( 0.5 );
// // Chan and Vese external term
ChanAndVeseExternalTermType::Pointer cvExternalTerm = ChanAndVeseExternalTermType::New();
cvExternalTerm->SetInput( input );
// Create term container (equation rhs)
TermContainerType::Pointer termContainer = TermContainerType::New();
termContainer->SetLevelSetContainer( levelSetContainer );
termContainer->SetInput( input );
termContainer->AddTerm( 0, cvInternalTerm );
termContainer->AddTerm( 1, cvExternalTerm );
// Create equation container
EquationContainerType::Pointer equationContainer = EquationContainerType::New();
equationContainer->SetLevelSetContainer( levelSetContainer );
equationContainer->AddEquation( 0, termContainer );
// Create stopping criteria
StoppingCriterionType::Pointer criterion = StoppingCriterionType::New();
criterion->SetNumberOfIterations( numberOfIterations );
// Create the visualizer
VisualizationType::Pointer visualizer = VisualizationType::New();
visualizer->SetInputImage( input );
visualizer->SetLevelSet( levelSet );
visualizer->SetScreenCapture( true );
// Create evolution class
LevelSetEvolutionType::Pointer evolution = LevelSetEvolutionType::New();
evolution->SetEquationContainer( equationContainer );
evolution->SetStoppingCriterion( criterion );
evolution->SetLevelSetContainer( levelSetContainer );
IterationUpdateCommandType::Pointer iterationUpdateCommand = IterationUpdateCommandType::New();
iterationUpdateCommand->SetFilterToUpdate( visualizer );
iterationUpdateCommand->SetUpdatePeriod( 5 );
evolution->AddObserver( itk::IterationEvent(), iterationUpdateCommand );
evolution->Update();
return EXIT_SUCCESS;
}