]> Creatis software - FrontAlgorithms.git/commitdiff
First commit
authorLeonardo Florez-Valencia <florez-l@javeriana.edu.co>
Tue, 9 Dec 2014 08:58:31 +0000 (09:58 +0100)
committerLeonardo Florez-Valencia <florez-l@javeriana.edu.co>
Tue, 9 Dec 2014 08:58:31 +0000 (09:58 +0100)
47 files changed:
CMakeLists.txt [new file with mode: 0644]
COMPILATION [new file with mode: 0644]
README [new file with mode: 0644]
appli/CMakeLists.txt [new file with mode: 0644]
appli/examples/CMakeLists.txt [new file with mode: 0644]
appli/examples/example_ImageAlgorithmDijkstra_00.cxx [new file with mode: 0644]
appli/examples/example_ImageAlgorithmDijkstra_01.cxx [new file with mode: 0644]
appli/examples/example_ImageAlgorithmDijkstra_02.cxx [new file with mode: 0644]
appli/examples/example_ImageAlgorithmFastMarching_00.cxx [new file with mode: 0644]
appli/examples/example_ImageAlgorithmFastMarching_01.cxx [new file with mode: 0644]
appli/examples/example_ImageAlgorithmRegionGrow_00.cxx [new file with mode: 0644]
appli/examples/example_ImageAlgorithmRegionGrow_01.cxx [new file with mode: 0644]
appli/examples/example_ImageAlgorithmRegionGrow_MultipleThresholds.cxx [new file with mode: 0644]
cmake/CMakeLists.txt [new file with mode: 0644]
cmake/FrontAlgorithmsConfig.cmake.in [new file with mode: 0644]
data/binary_test_2D_00.png [new file with mode: 0644]
data/ones_image.png [new file with mode: 0644]
lib/CMakeLists.txt [new file with mode: 0644]
lib/fpa/Base/Algorithm.h [new file with mode: 0644]
lib/fpa/Base/Algorithm.hxx [new file with mode: 0644]
lib/fpa/Base/Dijkstra.h [new file with mode: 0644]
lib/fpa/Base/Dijkstra.hxx [new file with mode: 0644]
lib/fpa/Base/Events.h [new file with mode: 0644]
lib/fpa/Base/FastMarching.h [new file with mode: 0644]
lib/fpa/Base/FastMarching.hxx [new file with mode: 0644]
lib/fpa/Base/Functors/InvertCostFunction.h [new file with mode: 0644]
lib/fpa/Base/RegionGrow.h [new file with mode: 0644]
lib/fpa/Base/RegionGrow.hxx [new file with mode: 0644]
lib/fpa/Base/TreeExtractor.h [new file with mode: 0644]
lib/fpa/Base/TreeExtractor.hxx [new file with mode: 0644]
lib/fpa/Common.cxx.in [new file with mode: 0644]
lib/fpa/Common.h [new file with mode: 0644]
lib/fpa/Image/Algorithm.h [new file with mode: 0644]
lib/fpa/Image/Algorithm.hxx [new file with mode: 0644]
lib/fpa/Image/Dijkstra.h [new file with mode: 0644]
lib/fpa/Image/FastMarching.h [new file with mode: 0644]
lib/fpa/Image/Functors/RegionGrowAllBelongsFunction.h [new file with mode: 0644]
lib/fpa/Image/Functors/RegionGrowThresholdFunction.h [new file with mode: 0644]
lib/fpa/Image/RegionGrow.h [new file with mode: 0644]
lib/fpa/Image/RegionGrowWithMultipleThresholds.h [new file with mode: 0644]
lib/fpa/Image/RegionGrowWithMultipleThresholds.hxx [new file with mode: 0644]
lib/fpa/VTK/Image2DObserver.h [new file with mode: 0644]
lib/fpa/VTK/Image2DObserver.hxx [new file with mode: 0644]
lib/fpa/VTK/Image3DObserver.h [new file with mode: 0644]
lib/fpa/VTK/Image3DObserver.hxx [new file with mode: 0644]
lib/fpa/VTK/ImageMPR.cxx [new file with mode: 0644]
lib/fpa/VTK/ImageMPR.h [new file with mode: 0644]

diff --git a/CMakeLists.txt b/CMakeLists.txt
new file mode 100644 (file)
index 0000000..cd777b0
--- /dev/null
@@ -0,0 +1,96 @@
+CMAKE_MINIMUM_REQUIRED(VERSION 2.6)
+
+# for CMake 2.6 corrected behaviour (see "cmake --help-policy CMP0003")
+IF(
+    COMMAND cmake_policy AND
+    ${CMAKE_MAJOR_VERSION} EQUAL 2 AND
+    ${CMAKE_MINOR_VERSION} GREATER 4
+    )
+  CMAKE_POLICY(SET CMP0003 NEW)
+  CMAKE_POLICY(SET CMP0005 NEW)
+  CMAKE_POLICY(SET CMP0011 NEW)
+  CMAKE_POLICY(SET CMP0012 NEW)
+ENDIF(
+  COMMAND cmake_policy AND
+  ${CMAKE_MAJOR_VERSION} EQUAL 2 AND
+  ${CMAKE_MINOR_VERSION} GREATER 4
+  )
+
+## ================
+## = Project name =
+## ================
+
+PROJECT(FrontAlgorithms)
+SET(FrontAlgorithms_MAJOR_VERSION   "0")
+SET(FrontAlgorithms_MINOR_VERSION   "0")
+SET(FrontAlgorithms_RELEASE_VERSION "1")
+SET(FrontAlgorithms_VERSION "${FrontAlgorithms_MAJOR_VERSION}.${FrontAlgorithms_MINOR_VERSION}.${FrontAlgorithms_RELEASE_VERSION}")
+
+## ===========
+## = Options =
+## ===========
+
+OPTION(BUILD_EXAMPLES "Build examples" OFF)
+OPTION(BUILD_SHARED_LIBS "Build shared libs" OFF)
+OPTION(USE_VTK "Build using VTK" OFF)
+
+IF(BUILD_SHARED_LIBS)
+  SET(LIB_TYPE SHARED)
+ELSE(BUILD_SHARED_LIBS)
+  SET(LIB_TYPE STATIC)
+ENDIF(BUILD_SHARED_LIBS)
+
+## ============
+## = Packages =
+## ============
+
+INCLUDE(GenerateExportHeader)
+
+FIND_PACKAGE(ITK REQUIRED)
+INCLUDE(${ITK_USE_FILE})
+
+IF(USE_VTK)
+  FIND_PACKAGE(VTK REQUIRED)
+  INCLUDE(${VTK_USE_FILE})
+ENDIF(USE_VTK)
+
+## ================================================
+## = Do not allow to build inside the source tree =
+## ================================================
+
+IF(PROJECT_BINARY_DIR STREQUAL ${PROJECT_SOURCE_DIR})
+  MESSAGE(FATAL_ERROR "Building in the source tree is not allowed")
+ENDIF(PROJECT_BINARY_DIR STREQUAL ${PROJECT_SOURCE_DIR})
+
+## =====================================
+## = Where to put executables and libs =
+## =====================================
+
+SET(EXECUTABLE_OUTPUT_PATH ${PROJECT_BINARY_DIR})
+SET(LIBRARY_OUTPUT_PATH ${PROJECT_BINARY_DIR})
+MARK_AS_ADVANCED(
+  CMAKE_BACKWARDS_COMPATIBILITY
+  EXECUTABLE_OUTPUT_PATH
+  LIBRARY_OUTPUT_PATH
+  )
+
+## ==============================
+## = Subdirs containing headers =
+## ==============================
+
+INCLUDE_DIRECTORIES(
+  ${PROJECT_SOURCE_DIR}/lib
+  ${PROJECT_BINARY_DIR}/lib
+  )
+
+## ===========================
+## = Subdirs containing code =
+## ===========================
+
+SUBDIRS(
+  cmake
+  lib
+  appli
+  )
+
+## eof - $RCSfile$
diff --git a/COMPILATION b/COMPILATION
new file mode 100644 (file)
index 0000000..1886dd8
--- /dev/null
@@ -0,0 +1,27 @@
+@description
+  The project uses CMake as compilation helper. It should compile on the three
+  major box flavors (linux, windows, mac). However, as of december 2014, it
+  has only been tested on linux Fedora 20 and Ubuntu 14.04.
+
+@prerequisites
+  1. CMake (>=2.8.12.2)
+
+  3. Insight Toolkit -ITK- (>=4.6.0)
+    3.1 Required cmake flags:
+      BUILD_SHARED_LIBS:BOOL=ON
+    3.2 [OPTIONAL] If USE_VTK=ON
+      Module_ITKVtkGlue:BOOL=ON
+
+  4. [OPTIONAL] Visualization Toolkit -VTK- (>=6.1.0)
+    4.1 Required cmake flags:
+      BUILD_SHARED_LIBS:BOOL=ON
+
+@cmake_flags
+  BUILD_EXAMPLES:BOOL
+    Build example applications? (most of them are command line)
+  CMAKE_BUILD_TYPE:STRING
+    Debug/Release?
+  CMAKE_INSTALL_PREFIX:STRING
+    Where to put installation products? (in windows this option has no use)
+
+## eof - $RCSfile$
diff --git a/README b/README
new file mode 100644 (file)
index 0000000..4ebdcfb
--- /dev/null
+++ b/README
@@ -0,0 +1,21 @@
+
+@project
+  FrontAlgorithms: Generic implementation of front propagation algorithms
+                   with some extra features
+
+@version
+  0.0.1 (2014-12-31)
+
+@authors
+  Maciej ORKISZ (maciej.orkisz@creatis.insa-lyon.fr)
+  Leonardo FLOREZ-VALENCIA (florez-l@javeriana.edu.co)
+  <Who else?>
+
+@description
+  <Write a fancy description and some historical notes>
+
+@license
+   <What do we put here?>
+
+
+## eof - $RCSfile$
diff --git a/appli/CMakeLists.txt b/appli/CMakeLists.txt
new file mode 100644 (file)
index 0000000..71b459c
--- /dev/null
@@ -0,0 +1,6 @@
+
+SUBDIRS(
+  examples
+  )
+
+## eof - $RCSfile$
diff --git a/appli/examples/CMakeLists.txt b/appli/examples/CMakeLists.txt
new file mode 100644 (file)
index 0000000..dff8fe0
--- /dev/null
@@ -0,0 +1,31 @@
+IF(BUILD_EXAMPLES)
+  SET(
+    APPLIS
+    example_ImageAlgorithmRegionGrow_00
+    example_ImageAlgorithmDijkstra_00
+    example_ImageAlgorithmFastMarching_00
+    )
+
+  FOREACH(APP ${APPLIS})
+    ADD_EXECUTABLE(${APP} ${APP}.cxx)
+    TARGET_LINK_LIBRARIES(${APP} FrontAlgorithms)
+  ENDFOREACH(APP)
+
+  IF(USE_VTK)
+    SET(
+      vtk_APPLIS
+      example_ImageAlgorithmRegionGrow_01
+      example_ImageAlgorithmRegionGrow_MultipleThresholds
+      example_ImageAlgorithmDijkstra_01
+      example_ImageAlgorithmDijkstra_02
+      example_ImageAlgorithmFastMarching_01
+      )
+
+    FOREACH(APP ${vtk_APPLIS})
+      ADD_EXECUTABLE(${APP} ${APP}.cxx)
+      TARGET_LINK_LIBRARIES(${APP} FrontAlgorithms)
+    ENDFOREACH(APP)
+  ENDIF(USE_VTK)
+ENDIF(BUILD_EXAMPLES)
+
+## eof - $RCSfile$
diff --git a/appli/examples/example_ImageAlgorithmDijkstra_00.cxx b/appli/examples/example_ImageAlgorithmDijkstra_00.cxx
new file mode 100644 (file)
index 0000000..86cff08
--- /dev/null
@@ -0,0 +1,44 @@
+#include <iostream>
+#include <itkImage.h>
+
+#include <fpa/Image/Dijkstra.h>
+
+// -------------------------------------------------------------------------
+const unsigned int Dim = 2;
+typedef unsigned char TPixel;
+typedef double        TCost;
+typedef itk::Image< TPixel, Dim > TImage;
+
+typedef fpa::Image::Dijkstra< TImage, TCost > TFrontAlgorithm;
+
+// -------------------------------------------------------------------------
+int main( int argc, char* argv[] )
+{
+  // Create a dummy image
+  TImage::SizeType imageSize;
+  imageSize.Fill( 100 );
+
+  TImage::SpacingType imageSpacing;
+  imageSpacing.Fill( 1 );
+
+  TImage::Pointer image = TImage::New( );
+  image->SetRegions( imageSize );
+  image->SetSpacing( imageSpacing );
+  image->Allocate( );
+  image->FillBuffer( TPixel( 1 ) );
+
+  // Seed
+  TImage::IndexType seed;
+  seed.Fill( 50 );
+
+  // Configure algorithm
+  TFrontAlgorithm::Pointer algorithm = TFrontAlgorithm::New( );
+  algorithm->AddSeed( seed, TCost( 0 ) );
+  algorithm->SetInput( image );
+  algorithm->SetNeighborhoodOrder( 1 );
+  algorithm->Update( );
+
+  return( 0 );
+}
+
+// eof - $RCSfile$
diff --git a/appli/examples/example_ImageAlgorithmDijkstra_01.cxx b/appli/examples/example_ImageAlgorithmDijkstra_01.cxx
new file mode 100644 (file)
index 0000000..b2614ed
--- /dev/null
@@ -0,0 +1,152 @@
+#include <iostream>
+#include <limits>
+#include <string>
+
+#include <itkImage.h>
+#include <itkImageFileReader.h>
+#include <itkImageToVTKImageFilter.h>
+
+#include <vtkImageActor.h>
+#include <vtkInteractorStyleImage.h>
+#include <vtkPointHandleRepresentation3D.h>
+#include <vtkProperty.h>
+#include <vtkRenderer.h>
+#include <vtkRenderWindow.h>
+#include <vtkRenderWindowInteractor.h>
+#include <vtkSeedRepresentation.h>
+#include <vtkSeedWidget.h>
+#include <vtkSmartPointer.h>
+
+#include <fpa/Image/Dijkstra.h>
+#include <fpa/VTK/Image2DObserver.h>
+
+// -------------------------------------------------------------------------
+const unsigned int Dim = 2;
+typedef unsigned char TPixel;
+typedef double TScalar;
+typedef itk::Image< TPixel, Dim >  TImage;
+typedef itk::ImageToVTKImageFilter< TImage > TVTKImage;
+
+typedef itk::ImageFileReader< TImage > TImageReader;
+typedef fpa::Image::Dijkstra< TImage, TScalar > TFrontAlgorithm;
+
+typedef
+fpa::VTK::Image2DObserver< TFrontAlgorithm, vtkRenderWindow >
+TObserver;
+
+// -------------------------------------------------------------------------
+int main( int argc, char* argv[] )
+{
+  if( argc < 2 )
+  {
+    std::cerr
+      << "Usage: " << argv[ 0 ]
+      << " input_image [stop_at_one_front]" << std::endl;
+    return( 1 );
+
+  } // fi
+  std::string input_image_fn = argv[ 1 ];
+  bool stop_at_one_front = false;
+  if( 2 < argc )
+    stop_at_one_front = ( std::atoi( argv[ 2 ] ) == 1 );
+
+  // Read image
+  TImageReader::Pointer input_image_reader = TImageReader::New( );
+  input_image_reader->SetFileName( input_image_fn );
+  try
+  {
+    input_image_reader->Update( );
+  }
+  catch( itk::ExceptionObject& err )
+  {
+    std::cerr << "Error caught: " << err << std::endl;
+    return( 1 );
+
+  } // yrt
+  TImage::ConstPointer input_image = input_image_reader->GetOutput( );
+
+  TVTKImage::Pointer vtk_image = TVTKImage::New( );
+  vtk_image->SetInput( input_image );
+  vtk_image->Update( );
+
+  // VTK visualization
+  vtkSmartPointer< vtkImageActor > actor =
+    vtkSmartPointer< vtkImageActor >::New( );
+  actor->SetInputData( vtk_image->GetOutput( ) );
+
+  vtkSmartPointer< vtkRenderer > renderer =
+    vtkSmartPointer< vtkRenderer >::New( );
+  renderer->SetBackground( 0.1, 0.2, 0.7 );
+  renderer->AddActor( actor );
+  vtkSmartPointer< vtkRenderWindow > window =
+    vtkSmartPointer< vtkRenderWindow >::New( );
+  window->SetSize( 800, 800 );
+  window->AddRenderer( renderer );
+
+  // VTK interaction
+  vtkSmartPointer< vtkInteractorStyleImage > imageStyle =
+    vtkSmartPointer< vtkInteractorStyleImage >::New( );
+  vtkSmartPointer< vtkRenderWindowInteractor > interactor =
+    vtkSmartPointer< vtkRenderWindowInteractor >::New( );
+  interactor->SetInteractorStyle( imageStyle );
+  window->SetInteractor( interactor );
+  window->Render( );
+
+  // Create the widget and its representation
+  vtkSmartPointer< vtkPointHandleRepresentation3D > handle =
+    vtkSmartPointer< vtkPointHandleRepresentation3D >::New( );
+  handle->GetProperty( )->SetColor( 1, 0, 0 );
+  vtkSmartPointer< vtkSeedRepresentation > rep =
+    vtkSmartPointer< vtkSeedRepresentation >::New( );
+  rep->SetHandleRepresentation( handle );
+
+  vtkSmartPointer< vtkSeedWidget > widget =
+    vtkSmartPointer< vtkSeedWidget >::New( );
+  widget->SetInteractor( interactor );
+  widget->SetRepresentation( rep );
+
+  // Let some interaction
+  interactor->Initialize( );
+  window->Render( );
+  widget->On( );
+  interactor->Start( );
+
+  // Configure observer
+  TObserver::Pointer obs = TObserver::New( );
+  obs->SetImage( input_image, window );
+
+  // Configure algorithm
+  TFrontAlgorithm::Pointer algorithm = TFrontAlgorithm::New( );
+  for( unsigned int s = 0; s < rep->GetNumberOfSeeds( ); s++ )
+  {
+    double pos[ 3 ];
+    rep->GetSeedWorldPosition( s, pos );
+
+    TImage::PointType pnt;
+    pnt[ 0 ] = pos[ 0 ];
+    pnt[ 1 ] = pos[ 1 ];
+
+    TImage::IndexType idx;
+    if( input_image->TransformPhysicalPointToIndex( pnt, idx ) )
+    {
+      algorithm->AddSeed( idx, 0 );
+      std::cout << " Seed --> " << idx << std::endl;
+
+    } // fi
+
+  } // rof
+  algorithm->AddObserver( itk::AnyEvent( ), obs );
+  algorithm->ThrowEventsOn( );
+  algorithm->SetInput( input_image );
+  algorithm->SetNeighborhoodOrder( 1 );
+  algorithm->SetStopAtOneFront( stop_at_one_front );
+  algorithm->Update( );
+
+  // One last interaction
+  window->Render( );
+  interactor->Start( );
+
+  return( 0 );
+}
+
+// eof - $RCSfile$
diff --git a/appli/examples/example_ImageAlgorithmDijkstra_02.cxx b/appli/examples/example_ImageAlgorithmDijkstra_02.cxx
new file mode 100644 (file)
index 0000000..728313e
--- /dev/null
@@ -0,0 +1,205 @@
+#include <iostream>
+#include <limits>
+#include <string>
+
+#include <itkImage.h>
+#include <itkImageFileReader.h>
+#include <itkRGBAPixel.h>
+#include <itkSignedDanielssonDistanceMapImageFilter.h>
+#include <itkImageToVTKImageFilter.h>
+
+#include <vtkImageActor.h>
+#include <vtkInteractorStyleImage.h>
+#include <vtkPointHandleRepresentation3D.h>
+#include <vtkProperty.h>
+#include <vtkRenderer.h>
+#include <vtkRenderWindow.h>
+#include <vtkRenderWindowInteractor.h>
+#include <vtkSeedRepresentation.h>
+#include <vtkSeedWidget.h>
+#include <vtkSmartPointer.h>
+
+#include <fpa/Image/Dijkstra.h>
+#include <fpa/Base/Functors/InvertCostFunction.h>
+#include <fpa/Base/TreeExtractor.h>
+#include <fpa/VTK/Image2DObserver.h>
+
+// -------------------------------------------------------------------------
+const unsigned int Dim = 2;
+typedef unsigned char TPixel;
+typedef itk::RGBAPixel< TPixel > TRGBAPixel;
+typedef double TScalar;
+typedef itk::Image< TPixel, Dim >  TImage;
+typedef itk::Image< TScalar, Dim > TDistanceMap;
+typedef itk::Image< TRGBAPixel, Dim > TRGBAImage;
+typedef itk::ImageToVTKImageFilter< TImage > TVTKImage;
+typedef itk::ImageToVTKImageFilter< TRGBAImage > TVTKRGBAImage;
+
+typedef itk::ImageFileReader< TImage > TImageReader;
+typedef
+itk::SignedDanielssonDistanceMapImageFilter< TImage, TDistanceMap >
+TDistanceFilter;
+
+typedef fpa::Image::Dijkstra< TDistanceMap, TScalar > TFrontAlgorithm;
+typedef fpa::Base::TreeExtractor< TFrontAlgorithm >   TExtractor;
+typedef fpa::Base::Functors::InvertCostFunction< TScalar > TScalarFunction;
+
+typedef
+fpa::VTK::Image2DObserver< TExtractor, vtkRenderWindow >
+TObserver;
+
+// -------------------------------------------------------------------------
+int main( int argc, char* argv[] )
+{
+  if( argc < 2 )
+  {
+    std::cerr
+      << "Usage: " << argv[ 0 ]
+      << " input_image [stop_at_one_front]" << std::endl;
+    return( 1 );
+
+  } // fi
+  std::string input_image_fn = argv[ 1 ];
+  bool stop_at_one_front = false;
+  if( 2 < argc )
+    stop_at_one_front = ( std::atoi( argv[ 2 ] ) == 1 );
+
+  // Read image
+  TImageReader::Pointer input_image_reader = TImageReader::New( );
+  input_image_reader->SetFileName( input_image_fn );
+  try
+  {
+    input_image_reader->Update( );
+  }
+  catch( itk::ExceptionObject& err )
+  {
+    std::cerr << "Error caught: " << err << std::endl;
+    return( 1 );
+
+  } // yrt
+    
+  TDistanceFilter::Pointer dist_filter = TDistanceFilter::New( );
+  dist_filter->SetInput( input_image_reader->GetOutput( ) );
+  dist_filter->InsideIsPositiveOn( );
+  dist_filter->SquaredDistanceOff( );
+  dist_filter->Update( );
+
+  TImage::ConstPointer input_image = input_image_reader->GetOutput( );
+
+  TVTKImage::Pointer vtk_image = TVTKImage::New( );
+  vtk_image->SetInput( input_image );
+  vtk_image->Update( );
+
+  // VTK visualization
+  vtkSmartPointer< vtkImageActor > actor =
+    vtkSmartPointer< vtkImageActor >::New( );
+  actor->SetInputData( vtk_image->GetOutput( ) );
+
+  vtkSmartPointer< vtkRenderer > renderer =
+    vtkSmartPointer< vtkRenderer >::New( );
+  renderer->SetBackground( 0.1, 0.2, 0.7 );
+  renderer->AddActor( actor );
+  vtkSmartPointer< vtkRenderWindow > window =
+    vtkSmartPointer< vtkRenderWindow >::New( );
+  window->SetSize( 800, 800 );
+  window->AddRenderer( renderer );
+
+  // VTK interaction
+  vtkSmartPointer< vtkInteractorStyleImage > imageStyle =
+    vtkSmartPointer< vtkInteractorStyleImage >::New( );
+  vtkSmartPointer< vtkRenderWindowInteractor > interactor =
+    vtkSmartPointer< vtkRenderWindowInteractor >::New( );
+  interactor->SetInteractorStyle( imageStyle );
+  window->SetInteractor( interactor );
+  window->Render( );
+
+  // Create the widget and its representation
+  vtkSmartPointer< vtkPointHandleRepresentation3D > handle =
+    vtkSmartPointer< vtkPointHandleRepresentation3D >::New( );
+  handle->GetProperty( )->SetColor( 1, 0, 0 );
+  vtkSmartPointer< vtkSeedRepresentation > rep =
+    vtkSmartPointer< vtkSeedRepresentation >::New( );
+  rep->SetHandleRepresentation( handle );
+
+  vtkSmartPointer< vtkSeedWidget > widget =
+    vtkSmartPointer< vtkSeedWidget >::New( );
+  widget->SetInteractor( interactor );
+  widget->SetRepresentation( rep );
+
+  // Let some interaction
+  interactor->Initialize( );
+  window->Render( );
+  widget->On( );
+  interactor->Start( );
+
+  // Configure observer
+  TObserver::Pointer obs = TObserver::New( );
+  obs->SetImage( dist_filter->GetOutput( ), window );
+
+  // Configure membership function
+  TScalarFunction::Pointer cost_function = TScalarFunction::New( );
+
+  // Configure algorithm
+  TExtractor::Pointer algorithm = TExtractor::New( );
+  for( unsigned int s = 0; s < rep->GetNumberOfSeeds( ); s++ )
+  {
+    double pos[ 3 ];
+    rep->GetSeedWorldPosition( s, pos );
+
+    TImage::PointType pnt;
+    pnt[ 0 ] = pos[ 0 ];
+    pnt[ 1 ] = pos[ 1 ];
+
+    TImage::IndexType idx;
+    if( input_image->TransformPhysicalPointToIndex( pnt, idx ) )
+    {
+      algorithm->AddSeed( idx, 0 );
+      if( s == 0 )
+        algorithm->SetRoot( idx );
+      else
+        algorithm->AddLeaf( idx );
+
+    } // fi
+
+  } // rof
+  algorithm->AddObserver( itk::AnyEvent( ), obs );
+  algorithm->ThrowEventsOn( );
+  algorithm->SetInput( dist_filter->GetOutput( ) );
+  algorithm->SetNeighborhoodOrder( 1 );
+  algorithm->SetStopAtOneFront( stop_at_one_front );
+  algorithm->SetCostConversion( cost_function );
+  algorithm->Update( );
+
+  // Create image with tree
+  TPixel transparent_color[] = { 0, 0, 0, 0 };
+  TPixel solid_color[] = { 255, 255, 0, 255 };
+  TRGBAImage::Pointer tree_image = TRGBAImage::New( );
+  tree_image->CopyInformation( input_image );
+  tree_image->SetRequestedRegionToLargestPossibleRegion( );
+  tree_image->SetBufferedRegion( tree_image->GetRequestedRegion( ) );
+  tree_image->Allocate( );
+  tree_image->FillBuffer( TRGBAPixel( transparent_color ) );
+
+  const TExtractor::TTree& tree = algorithm->GetTree( );
+  TExtractor::TTree::const_iterator tIt = tree.begin( );
+  for( ; tIt != tree.end( ); ++tIt )
+    tree_image->SetPixel( tIt->first, TRGBAPixel( solid_color ) );
+
+  TVTKRGBAImage::Pointer vtk_tree_image = TVTKRGBAImage::New( );
+  vtk_tree_image->SetInput( tree_image );
+  vtk_tree_image->Update( );
+
+  vtkSmartPointer< vtkImageActor > tree_actor =
+    vtkSmartPointer< vtkImageActor >::New( );
+  tree_actor->SetInputData( vtk_tree_image->GetOutput( ) );
+
+  renderer->AddActor( tree_actor );
+  window->Render( );
+
+  // One last interaction
+  interactor->Start( );
+
+  return( 0 );
+}
+
+// eof - $RCSfile$
diff --git a/appli/examples/example_ImageAlgorithmFastMarching_00.cxx b/appli/examples/example_ImageAlgorithmFastMarching_00.cxx
new file mode 100644 (file)
index 0000000..1573ce6
--- /dev/null
@@ -0,0 +1,44 @@
+#include <iostream>
+#include <itkImage.h>
+
+#include <fpa/Image/FastMarching.h>
+
+// -------------------------------------------------------------------------
+const unsigned int Dim = 2;
+typedef unsigned char TPixel;
+typedef double        TCost;
+typedef itk::Image< TPixel, Dim > TImage;
+
+typedef fpa::Image::FastMarching< TImage, TCost > TFrontAlgorithm;
+
+// -------------------------------------------------------------------------
+int main( int argc, char* argv[] )
+{
+  // Create a dummy image
+  TImage::SizeType imageSize;
+  imageSize.Fill( 100 );
+
+  TImage::SpacingType imageSpacing;
+  imageSpacing.Fill( 1 );
+
+  TImage::Pointer image = TImage::New( );
+  image->SetRegions( imageSize );
+  image->SetSpacing( imageSpacing );
+  image->Allocate( );
+  image->FillBuffer( TPixel( 1 ) );
+
+  // Seed
+  TImage::IndexType seed;
+  seed.Fill( 50 );
+
+  // Configure algorithm
+  TFrontAlgorithm::Pointer algorithm = TFrontAlgorithm::New( );
+  algorithm->AddSeed( seed, TCost( 0 ) );
+  algorithm->SetInput( image );
+  algorithm->SetNeighborhoodOrder( 1 );
+  algorithm->Update( );
+
+  return( 0 );
+}
+
+// eof - $RCSfile$
diff --git a/appli/examples/example_ImageAlgorithmFastMarching_01.cxx b/appli/examples/example_ImageAlgorithmFastMarching_01.cxx
new file mode 100644 (file)
index 0000000..cfa306f
--- /dev/null
@@ -0,0 +1,152 @@
+#include <iostream>
+#include <limits>
+#include <string>
+
+#include <itkImage.h>
+#include <itkImageFileReader.h>
+#include <itkImageToVTKImageFilter.h>
+
+#include <vtkImageActor.h>
+#include <vtkInteractorStyleImage.h>
+#include <vtkPointHandleRepresentation3D.h>
+#include <vtkProperty.h>
+#include <vtkRenderer.h>
+#include <vtkRenderWindow.h>
+#include <vtkRenderWindowInteractor.h>
+#include <vtkSeedRepresentation.h>
+#include <vtkSeedWidget.h>
+#include <vtkSmartPointer.h>
+
+#include <fpa/Image/FastMarching.h>
+#include <fpa/VTK/Image2DObserver.h>
+
+// -------------------------------------------------------------------------
+const unsigned int Dim = 2;
+typedef unsigned char TPixel;
+typedef double TScalar;
+typedef itk::Image< TPixel, Dim >  TImage;
+typedef itk::ImageToVTKImageFilter< TImage > TVTKImage;
+
+typedef itk::ImageFileReader< TImage > TImageReader;
+typedef fpa::Image::FastMarching< TImage, TScalar > TFrontAlgorithm;
+
+typedef
+fpa::VTK::Image2DObserver< TFrontAlgorithm, vtkRenderWindow >
+TObserver;
+
+// -------------------------------------------------------------------------
+int main( int argc, char* argv[] )
+{
+  if( argc < 2 )
+  {
+    std::cerr
+      << "Usage: " << argv[ 0 ]
+      << " input_image [stop_at_one_front]" << std::endl;
+    return( 1 );
+
+  } // fi
+  std::string input_image_fn = argv[ 1 ];
+  bool stop_at_one_front = false;
+  if( 2 < argc )
+    stop_at_one_front = ( std::atoi( argv[ 2 ] ) == 1 );
+
+  // Read image
+  TImageReader::Pointer input_image_reader = TImageReader::New( );
+  input_image_reader->SetFileName( input_image_fn );
+  try
+  {
+    input_image_reader->Update( );
+  }
+  catch( itk::ExceptionObject& err )
+  {
+    std::cerr << "Error caught: " << err << std::endl;
+    return( 1 );
+
+  } // yrt
+  TImage::ConstPointer input_image = input_image_reader->GetOutput( );
+
+  TVTKImage::Pointer vtk_image = TVTKImage::New( );
+  vtk_image->SetInput( input_image );
+  vtk_image->Update( );
+
+  // VTK visualization
+  vtkSmartPointer< vtkImageActor > actor =
+    vtkSmartPointer< vtkImageActor >::New( );
+  actor->SetInputData( vtk_image->GetOutput( ) );
+
+  vtkSmartPointer< vtkRenderer > renderer =
+    vtkSmartPointer< vtkRenderer >::New( );
+  renderer->SetBackground( 0.1, 0.2, 0.7 );
+  renderer->AddActor( actor );
+  vtkSmartPointer< vtkRenderWindow > window =
+    vtkSmartPointer< vtkRenderWindow >::New( );
+  window->SetSize( 800, 800 );
+  window->AddRenderer( renderer );
+
+  // VTK interaction
+  vtkSmartPointer< vtkInteractorStyleImage > imageStyle =
+    vtkSmartPointer< vtkInteractorStyleImage >::New( );
+  vtkSmartPointer< vtkRenderWindowInteractor > interactor =
+    vtkSmartPointer< vtkRenderWindowInteractor >::New( );
+  interactor->SetInteractorStyle( imageStyle );
+  window->SetInteractor( interactor );
+  window->Render( );
+
+  // Create the widget and its representation
+  vtkSmartPointer< vtkPointHandleRepresentation3D > handle =
+    vtkSmartPointer< vtkPointHandleRepresentation3D >::New( );
+  handle->GetProperty( )->SetColor( 1, 0, 0 );
+  vtkSmartPointer< vtkSeedRepresentation > rep =
+    vtkSmartPointer< vtkSeedRepresentation >::New( );
+  rep->SetHandleRepresentation( handle );
+
+  vtkSmartPointer< vtkSeedWidget > widget =
+    vtkSmartPointer< vtkSeedWidget >::New( );
+  widget->SetInteractor( interactor );
+  widget->SetRepresentation( rep );
+
+  // Let some interaction
+  interactor->Initialize( );
+  window->Render( );
+  widget->On( );
+  interactor->Start( );
+
+  // Configure observer
+  TObserver::Pointer obs = TObserver::New( );
+  obs->SetImage( input_image, window );
+
+  // Configure algorithm
+  TFrontAlgorithm::Pointer algorithm = TFrontAlgorithm::New( );
+  for( unsigned int s = 0; s < rep->GetNumberOfSeeds( ); s++ )
+  {
+    double pos[ 3 ];
+    rep->GetSeedWorldPosition( s, pos );
+
+    TImage::PointType pnt;
+    pnt[ 0 ] = pos[ 0 ];
+    pnt[ 1 ] = pos[ 1 ];
+
+    TImage::IndexType idx;
+    if( input_image->TransformPhysicalPointToIndex( pnt, idx ) )
+    {
+      algorithm->AddSeed( idx, 0 );
+      std::cout << " Seed --> " << idx << std::endl;
+
+    } // fi
+
+  } // rof
+  algorithm->AddObserver( itk::AnyEvent( ), obs );
+  algorithm->ThrowEventsOn( );
+  algorithm->SetInput( input_image );
+  algorithm->SetNeighborhoodOrder( 1 );
+  algorithm->SetStopAtOneFront( stop_at_one_front );
+  algorithm->Update( );
+
+  // One last interaction
+  window->Render( );
+  interactor->Start( );
+
+  return( 0 );
+}
+
+// eof - $RCSfile$
diff --git a/appli/examples/example_ImageAlgorithmRegionGrow_00.cxx b/appli/examples/example_ImageAlgorithmRegionGrow_00.cxx
new file mode 100644 (file)
index 0000000..149977a
--- /dev/null
@@ -0,0 +1,51 @@
+#include <iostream>
+#include <itkImage.h>
+
+#include <fpa/Image/RegionGrow.h>
+#include <fpa/Image/Functors/RegionGrowAllBelongsFunction.h>
+
+// -------------------------------------------------------------------------
+const unsigned int Dim = 2;
+typedef unsigned char TPixel;
+typedef itk::Image< TPixel, Dim > TImage;
+
+typedef fpa::Image::RegionGrow< TImage > TFrontAlgorithm;
+typedef
+fpa::Image::Functors::RegionGrowAllBelongsFunction< TImage >
+TMembershipFunction;
+
+// -------------------------------------------------------------------------
+int main( int argc, char* argv[] )
+{
+  // Create a dummy image
+  TImage::SizeType imageSize;
+  imageSize.Fill( 10 );
+
+  TImage::SpacingType imageSpacing;
+  imageSpacing.Fill( 1 );
+
+  TImage::Pointer image = TImage::New( );
+  image->SetRegions( imageSize );
+  image->SetSpacing( imageSpacing );
+  image->Allocate( );
+  image->FillBuffer( TPixel( 1 ) );
+
+  // Seed
+  TImage::IndexType seed;
+  seed.Fill( 5 );
+
+  // Configure membership function
+  TMembershipFunction::Pointer membership = TMembershipFunction::New( );
+
+  // Configure algorithm
+  TFrontAlgorithm::Pointer algorithm = TFrontAlgorithm::New( );
+  algorithm->AddSeed( seed, 0 );
+  algorithm->SetInput( image );
+  algorithm->SetNeighborhoodOrder( 1 );
+  algorithm->SetMembershipFunction( membership );
+  algorithm->Update( );
+
+  return( 0 );
+}
+
+// eof - $RCSfile$
diff --git a/appli/examples/example_ImageAlgorithmRegionGrow_01.cxx b/appli/examples/example_ImageAlgorithmRegionGrow_01.cxx
new file mode 100644 (file)
index 0000000..e41b65d
--- /dev/null
@@ -0,0 +1,157 @@
+#include <iostream>
+#include <limits>
+#include <string>
+
+#include <itkImage.h>
+#include <itkImageFileReader.h>
+#include <itkSignedDanielssonDistanceMapImageFilter.h>
+#include <itkImageToVTKImageFilter.h>
+
+#include <vtkImageActor.h>
+#include <vtkInteractorStyleImage.h>
+#include <vtkPointHandleRepresentation3D.h>
+#include <vtkProperty.h>
+#include <vtkRenderer.h>
+#include <vtkRenderWindow.h>
+#include <vtkRenderWindowInteractor.h>
+#include <vtkSeedRepresentation.h>
+#include <vtkSeedWidget.h>
+#include <vtkSmartPointer.h>
+
+#include <fpa/Image/RegionGrow.h>
+#include <fpa/Image/Functors/RegionGrowThresholdFunction.h>
+#include <fpa/VTK/Image2DObserver.h>
+
+// -------------------------------------------------------------------------
+const unsigned int Dim = 2;
+typedef unsigned char TPixel;
+typedef double TScalar;
+typedef itk::Image< TPixel, Dim >  TImage;
+typedef itk::ImageToVTKImageFilter< TImage > TVTKImage;
+
+typedef itk::ImageFileReader< TImage >   TImageReader;
+typedef fpa::Image::RegionGrow< TImage > TFrontAlgorithm;
+
+typedef
+fpa::Image::Functors::RegionGrowThresholdFunction< TImage >
+TMembershipFunction;
+
+typedef
+fpa::VTK::Image2DObserver< TFrontAlgorithm, vtkRenderWindow >
+TObserver;
+
+// -------------------------------------------------------------------------
+int main( int argc, char* argv[] )
+{
+  if( argc < 4 )
+  {
+    std::cerr
+      << "Usage: " << argv[ 0 ]
+      << " input_image lower_threshold upper_threshold" << std::endl;
+    return( 1 );
+
+  } // fi
+  std::string input_image_fn = argv[ 1 ];
+  TPixel lower_threshold = TPixel( std::atof( argv[ 2 ] ) );
+  TPixel upper_threshold = TPixel( std::atof( argv[ 3 ] ) );
+
+  // Read image
+  TImageReader::Pointer input_image_reader = TImageReader::New( );
+  input_image_reader->SetFileName( input_image_fn );
+  try
+  {
+    input_image_reader->Update( );
+  }
+  catch( itk::ExceptionObject& err )
+  {
+    std::cerr << "Error caught: " << err << std::endl;
+    return( 1 );
+
+  } // yrt
+  TImage::ConstPointer input_image = input_image_reader->GetOutput( );
+
+  TVTKImage::Pointer vtk_image = TVTKImage::New( );
+  vtk_image->SetInput( input_image );
+  vtk_image->Update( );
+
+  // VTK visualization
+  vtkSmartPointer< vtkImageActor > actor =
+    vtkSmartPointer< vtkImageActor >::New( );
+  actor->SetInputData( vtk_image->GetOutput( ) );
+
+  vtkSmartPointer< vtkRenderer > renderer =
+    vtkSmartPointer< vtkRenderer >::New( );
+  renderer->SetBackground( 0.1, 0.2, 0.7 );
+  renderer->AddActor( actor );
+  vtkSmartPointer< vtkRenderWindow > window =
+    vtkSmartPointer< vtkRenderWindow >::New( );
+  window->SetSize( 800, 800 );
+  window->AddRenderer( renderer );
+
+  // VTK interaction
+  vtkSmartPointer< vtkInteractorStyleImage > imageStyle =
+    vtkSmartPointer< vtkInteractorStyleImage >::New( );
+  vtkSmartPointer< vtkRenderWindowInteractor > interactor =
+    vtkSmartPointer< vtkRenderWindowInteractor >::New( );
+  interactor->SetInteractorStyle( imageStyle );
+  window->SetInteractor( interactor );
+  window->Render( );
+
+  // Create the widget and its representation
+  vtkSmartPointer< vtkPointHandleRepresentation3D > handle =
+    vtkSmartPointer< vtkPointHandleRepresentation3D >::New( );
+  handle->GetProperty( )->SetColor( 1, 0, 0 );
+  vtkSmartPointer< vtkSeedRepresentation > rep =
+    vtkSmartPointer< vtkSeedRepresentation >::New( );
+  rep->SetHandleRepresentation( handle );
+
+  vtkSmartPointer< vtkSeedWidget > widget =
+    vtkSmartPointer< vtkSeedWidget >::New( );
+  widget->SetInteractor( interactor );
+  widget->SetRepresentation( rep );
+
+  // Let some interaction
+  interactor->Initialize( );
+  window->Render( );
+  widget->On( );
+  interactor->Start( );
+
+  // Configure observer
+  TObserver::Pointer obs = TObserver::New( );
+  obs->SetImage( input_image, window );
+
+  // Configure membership function
+  TMembershipFunction::Pointer membership = TMembershipFunction::New( );
+  membership->SetLowerThreshold( lower_threshold );
+  membership->SetUpperThreshold( upper_threshold );
+
+  // Configure algorithm
+  TFrontAlgorithm::Pointer algorithm = TFrontAlgorithm::New( );
+  for( unsigned int s = 0; s < rep->GetNumberOfSeeds( ); s++ )
+  {
+    double pos[ 3 ];
+    rep->GetSeedWorldPosition( s, pos );
+
+    TImage::PointType pnt;
+    pnt[ 0 ] = TScalar( pos[ 0 ] );
+    pnt[ 1 ] = TScalar( pos[ 1 ] );
+
+    TImage::IndexType idx;
+    if( input_image->TransformPhysicalPointToIndex( pnt, idx ) )
+      algorithm->AddSeed( idx, 0 );
+
+  } // rof
+  algorithm->AddObserver( itk::AnyEvent( ), obs );
+  algorithm->ThrowEventsOn( );
+  algorithm->SetInput( input_image );
+  algorithm->SetNeighborhoodOrder( 1 );
+  algorithm->SetMembershipFunction( membership );
+  algorithm->Update( );
+
+  // One last interaction
+  interactor->Start( );
+
+  return( 0 );
+}
+
+// eof - $RCSfile$
diff --git a/appli/examples/example_ImageAlgorithmRegionGrow_MultipleThresholds.cxx b/appli/examples/example_ImageAlgorithmRegionGrow_MultipleThresholds.cxx
new file mode 100644 (file)
index 0000000..165dafb
--- /dev/null
@@ -0,0 +1,118 @@
+#include <iostream>
+#include <limits>
+#include <set>
+#include <string>
+
+#include <itkImage.h>
+#include <itkImageFileReader.h>
+#include <itkImageToVTKImageFilter.h>
+
+#include <vtkImageActor.h>
+#include <vtkProperty.h>
+#include <vtkRenderer.h>
+#include <vtkRenderWindow.h>
+#include <vtkRenderWindowInteractor.h>
+#include <vtkSmartPointer.h>
+#include <vtkSphereSource.h>
+
+#include <fpa/Image/RegionGrowWithMultipleThresholds.h>
+#include <fpa/VTK/ImageMPR.h>
+#include <fpa/VTK/Image3DObserver.h>
+
+// -------------------------------------------------------------------------
+const unsigned int Dim = 3;
+typedef short TPixel;
+typedef double TScalar;
+typedef itk::Image< TPixel, Dim >  TImage;
+typedef itk::ImageToVTKImageFilter< TImage > TVTKImage;
+
+typedef itk::ImageFileReader< TImage >   TImageReader;
+
+typedef
+fpa::Image::RegionGrowWithMultipleThresholds< TImage >
+TFrontAlgorithm;
+typedef
+fpa::VTK::Image3DObserver< TFrontAlgorithm, vtkRenderWindow >
+TObserver;
+
+// -------------------------------------------------------------------------
+int main( int argc, char* argv[] )
+{
+  if( argc < 8 )
+  {
+    std::cerr
+      << "Usage: " << argv[ 0 ]
+      << " input_image s_x s_y s_z thr_0 thr_1 n_samples" << std::endl;
+    return( 1 );
+
+  } // fi
+  std::string input_image_fn = argv[ 1 ];
+  TImage::PointType seed_pnt;
+  seed_pnt[ 0 ] = std::atof( argv[ 2 ] );
+  seed_pnt[ 1 ] = std::atof( argv[ 3 ] );
+  seed_pnt[ 2 ] = std::atof( argv[ 4 ] );
+  TPixel thr_0 = TPixel( std::atof( argv[ 5 ] ) );
+  TPixel thr_1 = TPixel( std::atof( argv[ 6 ] ) );
+  unsigned int n_samples = std::atoi( argv[ 7 ] );
+  
+  // Read image
+  TImageReader::Pointer input_image_reader = TImageReader::New( );
+  input_image_reader->SetFileName( input_image_fn );
+  try
+  {
+    input_image_reader->Update( );
+  }
+  catch( itk::ExceptionObject& err )
+  {
+    std::cerr << "Error caught: " << err << std::endl;
+    return( 1 );
+
+  } // yrt
+  TImage::ConstPointer input_image = input_image_reader->GetOutput( );
+  TImage::SpacingType spac = input_image->GetSpacing( );
+  double min_spac = spac[ 0 ];
+  for( unsigned int d = 1; d < Dim; d++ )
+    min_spac = ( spac[ d ] < min_spac )? spac[ d ]: min_spac;
+
+  TVTKImage::Pointer vtk_image = TVTKImage::New( );
+  vtk_image->SetInput( input_image );
+  vtk_image->Update( );
+
+  vtkSmartPointer< vtkSphereSource > seed =
+    vtkSmartPointer< vtkSphereSource >::New( );
+  seed->SetCenter( seed_pnt[ 0 ], seed_pnt[ 1 ], seed_pnt[ 2 ] );
+  seed->SetRadius( min_spac * double( 5 ) );
+  seed->Update( );
+
+  fpa::VTK::ImageMPR view;
+  view.SetBackground( 0.3, 0.2, 0.8 );
+  view.SetSize( 800, 800 );
+  view.SetImage( vtk_image->GetOutput( ) );
+  view.AddPolyData( seed->GetOutput( ), 1, 0, 0 );
+  view.Start( );
+
+  // Configure observer
+  TObserver::Pointer obs = TObserver::New( );
+  obs->SetImage( input_image, view.GetWindow( ) );
+
+  // Configure algorithm
+  TFrontAlgorithm::Pointer algorithm = TFrontAlgorithm::New( );
+  algorithm->AddThresholds( thr_0, thr_1, n_samples );
+
+  TImage::IndexType seed_idx;
+  input_image->TransformPhysicalPointToIndex( seed_pnt, seed_idx );
+  algorithm->AddSeed( seed_idx, 0 );
+
+  algorithm->AddObserver( itk::AnyEvent( ), obs );
+  algorithm->ThrowEventsOn( );
+  algorithm->SetInput( input_image );
+  algorithm->SetNeighborhoodOrder( 1 );
+  algorithm->SetDerivativeThreshold( double( 3 ) );
+  algorithm->Update( );
+
+  view.Start( );
+
+  return( 0 );
+}
+
+// eof - $RCSfile$
diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt
new file mode 100644 (file)
index 0000000..0123559
--- /dev/null
@@ -0,0 +1,7 @@
+CONFIGURE_FILE(
+  FrontAlgorithmsConfig.cmake.in
+  ${PROJECT_BINARY_DIR}/FrontAlgorithmsConfig.cmake
+  @ONLY
+  )
+
+## eof - $RCSfile$
diff --git a/cmake/FrontAlgorithmsConfig.cmake.in b/cmake/FrontAlgorithmsConfig.cmake.in
new file mode 100644 (file)
index 0000000..da3d3e8
--- /dev/null
@@ -0,0 +1,43 @@
+FIND_PATH(
+  FrontAlgorithms_INCLUDE_DIR1
+  fpa/Base/Algorithm.h
+  PATHS
+  /usr/include
+  /usr/local/include
+  @PROJECT_SOURCE_DIR@/lib
+  @PROJECT_BINARY_DIR@/lib
+  @CMAKE_INSTALL_PREFIX@/include
+  )
+
+FIND_PATH(
+  FrontAlgorithms_INCLUDE_DIR2
+  fpa/FrontAlgorithms_Export.h
+  PATHS
+  /usr/include
+  /usr/local/include
+  @PROJECT_SOURCE_DIR@/lib
+  @PROJECT_BINARY_DIR@/lib
+  @CMAKE_INSTALL_PREFIX@/include
+  )
+
+INCLUDE_DIRECTORIES(
+  ${FrontAlgorithms_INCLUDE_DIR1}
+  ${FrontAlgorithms_INCLUDE_DIR2}
+  )
+
+FIND_LIBRARY(
+  FrontAlgorithms_LIBRARY_NAME
+  FrontAlgorithms
+  PATHS
+  /usr/lib
+  /usr/local/lib
+  @PROJECT_BINARY_DIR@
+  @CMAKE_INSTALL_PREFIX@/lib
+  )
+
+SET(
+  FrontAlgorithms_LIBRARIES
+  ${FrontAlgorithms_LIBRARY_NAME}
+  )
+
+## eof - $RCSfile$
diff --git a/data/binary_test_2D_00.png b/data/binary_test_2D_00.png
new file mode 100644 (file)
index 0000000..fbe9d39
Binary files /dev/null and b/data/binary_test_2D_00.png differ
diff --git a/data/ones_image.png b/data/ones_image.png
new file mode 100644 (file)
index 0000000..f7afaab
Binary files /dev/null and b/data/ones_image.png differ
diff --git a/lib/CMakeLists.txt b/lib/CMakeLists.txt
new file mode 100644 (file)
index 0000000..f1821ee
--- /dev/null
@@ -0,0 +1,67 @@
+SET(LIB_NAME FrontAlgorithms)
+
+## ================
+## = Source files =
+## ================
+
+CONFIGURE_FILE(
+  fpa/Common.cxx.in
+  ${PROJECT_BINARY_DIR}/lib/fpa/Common.cxx
+  )
+
+FILE(GLOB ${LIB_NAME}_HEADERS "fpa/*.h" "fpa/*.hxx")
+FILE(GLOB ${LIB_NAME}_BASE_HEADERS "fpa/Base/*.h" "fpa/Base/*.hxx")
+FILE(GLOB ${LIB_NAME}_IMAGE_HEADERS "fpa/Image/*.h" "fpa/Image/*.hxx")
+
+FILE(GLOB ${LIB_NAME}_SOURCES "fpa/*.cxx")
+FILE(GLOB ${LIB_NAME}_BASE_SOURCES "fpa/Base/*.cxx")
+FILE(GLOB ${LIB_NAME}_IMAGE_SOURCES "fpa/Image/*.cxx")
+
+IF(USE_VTK)
+  FILE(GLOB ${LIB_NAME}_VTK_HEADERS "fpa/VTK/*.h" "fpa/VTK/*.hxx")
+  FILE(GLOB ${LIB_NAME}_VTK_SOURCES "fpa/VTK/*.cxx")
+ENDIF(USE_VTK)
+
+SET(
+  ${LIB_NAME}_ALL_SOURCES
+  ${PROJECT_BINARY_DIR}/lib/fpa/Common.cxx
+  ${${LIB_NAME}_SOURCES}
+  ${${LIB_NAME}_BASE_SOURCES}
+  ${${LIB_NAME}_IMAGE_SOURCES}
+  ${${LIB_NAME}_VTK_SOURCES}
+  )
+
+## =============
+## = Libraries =
+## =============
+
+SET(
+  ${LIB_NAME}_LINK_LIBRARIES
+  ${ITK_LIBRARIES}
+  ${VTK_LIBRARIES}
+  )
+
+## =====================
+## = Compilation rules =
+## =====================
+
+ADD_LIBRARY(
+  ${LIB_NAME}
+  ${LIB_TYPE}
+  ${${LIB_NAME}_ALL_SOURCES}
+  )
+GENERATE_EXPORT_HEADER(
+  ${LIB_NAME}
+  BASE_NAME ${LIB_NAME}
+  EXPORT_MACRO_NAME ${LIB_NAME}_EXPORT
+  EXPORT_FILE_NAME ${PROJECT_BINARY_DIR}/lib/fpa/${LIB_NAME}_Export.h
+  STATIC_DEFINE ${LIB_NAME}_BUILT_AS_STATIC
+  )
+TARGET_LINK_LIBRARIES(
+  ${LIB_NAME}
+  ${${LIB_NAME}_LINK_LIBRARIES}
+  ${ITK_LIBRARIES}
+  ${VTK_LIBRARIES}
+  )
+
+## eof - $RCSfile$
diff --git a/lib/fpa/Base/Algorithm.h b/lib/fpa/Base/Algorithm.h
new file mode 100644 (file)
index 0000000..187a281
--- /dev/null
@@ -0,0 +1,153 @@
+#ifndef __FPA__BASE__ALGORITHM__H__
+#define __FPA__BASE__ALGORITHM__H__
+
+#include <map>
+#include <vector>
+#include <utility>
+#include <fpa/Base/Events.h>
+
+namespace fpa
+{
+  namespace Base
+  {
+    /**
+     * Base front propagation algorithm. From a series of start seeds with
+     * costs, a priority queue is filled and emptied updating costs. Each
+     * vertex could be marked as "visited", "in the front", "not yet there"
+     * or "freezed".
+     *
+     * @param T Traits used for this algorithm
+     */
+    template< class T, class B >
+    class Algorithm
+      : public B
+    {
+    public:
+      // Standard class typdedefs
+      typedef Algorithm                       Self;
+      typedef B                               Superclass;
+      typedef itk::SmartPointer< Self >       Pointer;
+      typedef itk::SmartPointer< const Self > ConstPointer;
+
+      /// Templated types
+      typedef T TTraits;
+      typedef B TBaseFilter;
+      typedef typename T::TCost        TCost;
+      typedef typename T::TResult      TResult;
+      typedef typename T::TVertex      TVertex;
+      typedef typename T::TVertexValue TVertexValue;
+
+    protected:
+      typedef typename T::TFrontId   _TFrontId;
+      typedef typename T::TNode      _TNode;
+      typedef typename T::TNodes     _TNodes;
+      typedef typename T::TVertexCmp _TVertexCmp;
+
+      typedef std::map< TVertex, _TNode, _TVertexCmp >   _TMarks;
+
+      typedef std::pair< TVertex, bool >         _TCollision;
+      typedef std::vector< _TCollision >         _TCollisionSitesRow;
+      typedef std::vector< _TCollisionSitesRow > _TCollisionSites;
+
+    public:
+      typedef BaseEvent< _TNode >      TEvent;
+      typedef FrontEvent< _TNode >     TFrontEvent;
+      typedef MarkEvent< _TNode >      TMarkEvent;
+      typedef CollisionEvent< _TNode > TCollisionEvent;
+      typedef EndEvent< _TNode >       TEndEvent;
+
+    public:
+      itkTypeMacro( Algorithm, itkProcessObject );
+
+      itkBooleanMacro( StopAtOneFront );
+      itkBooleanMacro( ThrowEvents );
+
+      itkGetConstMacro( StopAtOneFront, bool );
+      itkGetConstMacro( ThrowEvents, bool );
+
+      itkSetMacro( StopAtOneFront, bool );
+      itkSetMacro( ThrowEvents, bool );
+
+    public:
+      virtual void InvokeEvent( const itk::EventObject& e );
+      virtual void InvokeEvent( const itk::EventObject& e ) const;
+
+      /// Seeds manipulation
+      void AddSeed( const TVertex& s, const TResult& v );
+      const TVertex& GetSeed( const unsigned long& i ) const;
+      void ClearSeeds( );
+      unsigned long GetNumberOfSeeds( ) const;
+
+    protected:
+      Algorithm( );
+      virtual ~Algorithm( );
+
+      /// itk::ProcessObject
+      virtual void GenerateData( );
+
+      /// Base interface
+      virtual void _BeforeLoop         ( );
+      virtual void _AfterLoop          ( );
+      virtual void _Loop               ( );
+      virtual bool _CheckCollisions    ( const _TNode& a, const _TNode& b );
+      virtual bool _CheckStopCondition ( );
+      virtual bool _UpdateResult       ( _TNode& n );
+
+      /// Marks management
+      virtual      void _InitializeMarks ( );
+      virtual      bool _IsMarked        ( const TVertex& v ) const;
+      virtual _TFrontId _FrontId         ( const TVertex& v ) const;
+      virtual   TVertex _Parent          ( const TVertex& v ) const;
+      virtual      void _Mark            ( const _TNode& n );
+
+      /// Pure virtual interface: vertices
+      virtual unsigned long _NumberOfVertices ( ) const = 0;
+      virtual  TVertexValue _Value            ( const TVertex& v ) const = 0;
+      virtual       TResult _Result           ( const TVertex& v ) const = 0;
+
+      /// Pure virtual interface: edges
+      virtual double _Norm ( const TVertex& a, const TVertex& b ) const = 0;
+      virtual   bool _Edge ( const TVertex& a, const TVertex& b ) const = 0;
+      virtual  TCost _Cost ( const TVertex& a, const TVertex& b ) const = 0;
+
+      /// Pure virtual interface: neighborhood
+      virtual void _Neighs      ( const _TNode& n, _TNodes& N ) const = 0;
+      virtual bool _UpdateNeigh ( _TNode& nn, const _TNode& n ) = 0;
+      virtual void _NeighsInDim ( const _TNode& n,
+                                  const unsigned int& d,
+                                  _TNodes& N ) = 0;
+
+      /// Pure virtual interface: queue
+      virtual   void _InitializeQueue ( ) = 0;
+      virtual   bool _IsQueueEmpty    ( ) const = 0;
+      virtual   void _QueuePush       ( const _TNode& n ) = 0;
+      virtual _TNode _QueuePop        ( ) = 0;
+      virtual   void _QueueClear      ( ) = 0;
+
+      /// Pure virtual interface: results
+      virtual void _InitializeResults ( ) = 0;
+
+    private:
+      // Purposely not implemented
+      Algorithm( const Self& );
+      void operator=( const Self& );
+
+    protected:
+      bool m_StopAtOneFront;
+      bool m_ThrowEvents;
+
+      _TNodes m_Seeds;
+      _TMarks m_Marks;
+
+      _TCollisionSites m_CollisionSites;
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#include <fpa/Base/Algorithm.hxx>
+
+#endif // __FPA__BASE__ALGORITHM__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Base/Algorithm.hxx b/lib/fpa/Base/Algorithm.hxx
new file mode 100644 (file)
index 0000000..286343f
--- /dev/null
@@ -0,0 +1,300 @@
+#ifndef __FPA__BASE__ALGORITHM__HXX__
+#define __FPA__BASE__ALGORITHM__HXX__
+
+#include <algorithm>
+#include <limits>
+#include <queue>
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+void fpa::Base::Algorithm< T, B >::
+InvokeEvent( const itk::EventObject& e )
+{
+  if( this->m_ThrowEvents )
+    this->Superclass::InvokeEvent( e );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+void fpa::Base::Algorithm< T, B >::
+InvokeEvent( const itk::EventObject& e ) const
+{
+  if( this->m_ThrowEvents )
+    this->Superclass::InvokeEvent( e );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+void fpa::Base::Algorithm< T, B >::
+AddSeed( const TVertex& s, const TResult& v )
+{
+  this->m_Seeds.push_back( _TNode( s, v, this->m_Seeds.size( ) ) );
+  this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+const typename fpa::Base::Algorithm< T, B >::
+TVertex& fpa::Base::Algorithm< T, B >::
+GetSeed( const unsigned long& i ) const
+{
+  return( this->m_Seeds[ i ].Vertex );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+void fpa::Base::Algorithm< T, B >::
+ClearSeeds( )
+{
+  this->m_Seeds.clear( );
+  this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+unsigned long fpa::Base::Algorithm< T, B >::
+GetNumberOfSeeds( ) const
+{
+  return( this->m_Seeds.size( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+fpa::Base::Algorithm< T, B >::
+Algorithm( )
+  : Superclass( ),
+    m_StopAtOneFront( false ),
+    m_ThrowEvents( false )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+fpa::Base::Algorithm< T, B >::
+~Algorithm( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+void fpa::Base::Algorithm< T, B >::
+GenerateData( )
+{
+  this->AllocateOutputs( );
+
+  unsigned long N = this->m_Seeds.size( );
+  if( N == 0 )
+    return;
+  this->m_CollisionSites.clear( );
+  this->m_CollisionSites.
+    resize( N, _TCollisionSitesRow( N, _TCollision( TVertex( ), false ) ) );
+  this->_BeforeLoop( );
+  this->_Loop( );
+  this->_AfterLoop( );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+void fpa::Base::Algorithm< T, B >::
+_BeforeLoop( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+void fpa::Base::Algorithm< T, B >::
+_AfterLoop( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+void fpa::Base::Algorithm< T, B >::
+_Loop( )
+{
+  this->_InitializeMarks( );
+  this->_InitializeResults( );
+  this->_InitializeQueue( );
+  while( !( this->_IsQueueEmpty( ) ) )
+  {
+    _TNode n = this->_QueuePop( );
+    if( this->_IsMarked( n.Vertex ) )
+      continue;
+
+    this->_Mark( n );
+    this->InvokeEvent( TMarkEvent( n ) );
+    if( this->_UpdateResult( n ) )
+    {
+      if( !( this->_CheckStopCondition( ) ) )
+      {
+        _TNodes N;
+        this->_Neighs( n, N );
+        typename _TNodes::iterator nnIt = N.begin( );
+        while( nnIt != N.end( ) )
+        {
+          if( this->_IsMarked( nnIt->Vertex ) )
+          {
+            // Update real front identifier
+            nnIt->FrontId = this->_FrontId( nnIt->Vertex );
+
+            // Update collisions
+            if( this->_CheckCollisions( n, *nnIt ) )
+            {
+              if( this->m_StopAtOneFront )
+              {
+                this->_QueueClear( );
+                nnIt = N.end( );
+
+              } // fi
+
+            } // fi
+          }
+          else
+          {
+            if( this->_UpdateNeigh( *nnIt, n ) )
+            {
+              nnIt->Parent = n.Vertex;
+              this->_QueuePush( *nnIt );
+              this->InvokeEvent( TFrontEvent( *nnIt ) );
+
+            } // fi
+
+          } // fi
+          if( nnIt != N.end( ) )
+            nnIt++;
+
+        } // elihw
+      }
+      else
+        this->_QueueClear( );
+
+    } // fi
+
+  } // elihw
+  this->InvokeEvent( TEndEvent( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+bool fpa::Base::Algorithm< T, B >::
+_CheckCollisions( const _TNode& a, const _TNode& b )
+{
+  bool ret = false;
+  if( a.FrontId != b.FrontId )
+  {
+    // Mark collision, if it is new
+    bool exists = this->m_CollisionSites[ a.FrontId ][ b.FrontId ].second;
+    exists     &= this->m_CollisionSites[ b.FrontId ][ a.FrontId ].second;
+    if( !exists )
+    {
+      this->m_CollisionSites[ a.FrontId ][ b.FrontId ].first = a.Vertex;
+      this->m_CollisionSites[ a.FrontId ][ b.FrontId ].second = true;
+      this->m_CollisionSites[ b.FrontId ][ a.FrontId ].first = b.Vertex;
+      this->m_CollisionSites[ b.FrontId ][ a.FrontId ].second = true;
+
+      // Stop if one front is desired
+      if( this->m_StopAtOneFront )
+      {
+        // Perform a depth-first iteration on front graph
+        unsigned long N = this->GetNumberOfSeeds( );
+        unsigned long C = 0;
+        std::vector< bool > m( N, false );
+        std::queue< unsigned long > q;
+        q.push( 0 );
+        while( !q.empty( ) )
+        {
+          unsigned long f = q.front( );
+          q.pop( );
+
+          if( m[ f ] )
+            continue;
+          m[ f ] = true;
+          C++;
+
+          for( unsigned int n = 0; n < N; ++n )
+            if( this->m_CollisionSites[ f ][ n ].second && !m[ n ] )
+              q.push( n );
+
+        } // elihw
+        ret = ( C == N );
+
+      } // fi
+
+    } // fi
+
+  } // fi
+  return( ret );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+bool fpa::Base::Algorithm< T, B >::
+_CheckStopCondition( )
+{
+  return( false );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+bool fpa::Base::Algorithm< T, B >::
+_UpdateResult( _TNode& n )
+{
+  return( true );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+void fpa::Base::Algorithm< T, B >::
+_InitializeMarks( )
+{
+  this->m_Marks.clear( );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+bool  fpa::Base::Algorithm< T, B >::
+_IsMarked( const TVertex& v ) const
+{
+  return( this->m_Marks.find( v ) != this->m_Marks.end( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+typename fpa::Base::Algorithm< T, B >::
+_TFrontId fpa::Base::Algorithm< T, B >::
+_FrontId( const TVertex& v ) const
+{
+  typename _TMarks::const_iterator mIt = this->m_Marks.find( v );
+  if( mIt != this->m_Marks.end( ) )
+    return( mIt->second.FrontId );
+  else
+    return( std::numeric_limits< _TFrontId >::max( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+typename fpa::Base::Algorithm< T, B >::
+TVertex fpa::Base::Algorithm< T, B >::
+_Parent( const TVertex& v ) const
+{
+  typename _TMarks::const_iterator mIt = this->m_Marks.find( v );
+  if( mIt == this->m_Marks.end( ) )
+  {
+    TVertex v;
+    return( v );
+  }
+  else
+    return( mIt->second.Parent );
+}
+
+// -------------------------------------------------------------------------
+template< class T, class B >
+void  fpa::Base::Algorithm< T, B >::
+_Mark( const _TNode& n )
+{
+  this->m_Marks[ n.Vertex ] = n;
+}
+
+#endif // __FPA__BASE__ALGORITHM__HXX__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Base/Dijkstra.h b/lib/fpa/Base/Dijkstra.h
new file mode 100644 (file)
index 0000000..a6212cf
--- /dev/null
@@ -0,0 +1,120 @@
+#ifndef __FPA__BASE__DIJKSTRA__H__
+#define __FPA__BASE__DIJKSTRA__H__
+
+#include <vector>
+#include <fpa/Base/Algorithm.h>
+
+namespace fpa
+{
+  namespace Base
+  {
+    /**
+     */
+    template< class V, class C, class VV, class VC >
+    class DijkstraTraits
+    {
+    public:
+      typedef V TVertex;
+      typedef C TResult;
+      typedef C TCost;
+      typedef VV TVertexValue;
+      typedef VC TVertexCmp;
+      typedef long TFrontId;
+
+      class TNode
+      {
+      public:
+        TNode( )
+          : Cost( 0 )
+          { }
+        TNode( const TVertex& v, const TFrontId& f )
+          : Vertex( v ),
+            Parent( v ),
+            Result( TResult( 0 ) ),
+            FrontId( f ),
+            Cost( TCost( 0 ) )
+          { }
+        TNode( const TVertex& v, const TResult& r, const TFrontId& f )
+          : Vertex( v ),
+            Parent( v ),
+            Result( r ),
+            FrontId( f ),
+            Cost( TCost( 0 ) )
+          { }
+        virtual ~TNode( )
+          { }
+
+        // NOTE: stl::heaps work as maximum priority queues
+        bool operator<( const TNode& other ) const
+          { return( other.Cost < this->Cost ); }
+
+        TVertex Vertex;
+        TVertex Parent;
+        TResult Result;
+        TFrontId FrontId;
+        TCost Cost;
+      };
+
+      typedef std::vector< TNode > TNodes;
+    };
+
+    /**
+     * Dijkstra is a front propagation algorithm that minimizes costs
+     */
+    template< class V, class C, class VV, class VC, class B >
+    class Dijkstra
+      : public Algorithm< DijkstraTraits< V, C, VV, VC >, B >
+    {
+    public:
+      // Templated types
+      typedef V TVertex;
+      typedef C TCost;
+      typedef VV TVertexValue;
+      typedef B  TBaseFilter;
+      typedef DijkstraTraits< V, C, VV, VC > TTraits;
+
+      // Standard class typdedefs
+      typedef Dijkstra                        Self;
+      typedef Algorithm< TTraits, B >         Superclass;
+      typedef itk::SmartPointer< Self >       Pointer;
+      typedef itk::SmartPointer< const Self > ConstPointer;
+
+    protected:
+      typedef typename TTraits::TFrontId  _TFrontId;
+      typedef typename TTraits::TNode     _TNode;
+      typedef typename TTraits::TNodes    _TNodes;
+
+      typedef std::vector< _TNode > _TQueue;
+
+    public:
+      itkTypeMacro( Dijkstra, Base );
+
+    protected:
+      Dijkstra( );
+      virtual ~Dijkstra( );
+
+      virtual   void _InitializeQueue ( );
+      virtual   bool _IsQueueEmpty    ( ) const;
+      virtual   void _QueuePush       ( const _TNode& n );
+      virtual _TNode _QueuePop        ( );
+      virtual   void _QueueClear      ( );
+      virtual   bool _UpdateNeigh     ( _TNode& nn, const _TNode& n );
+
+    private:
+      // Purposely not implemented
+      Dijkstra( const Self& );
+      void operator=( const Self& );
+
+    private:
+      _TQueue m_Queue;
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#include <fpa/Base/Dijkstra.hxx>
+
+#endif // __FPA__BASE__DIJKSTRA__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Base/Dijkstra.hxx b/lib/fpa/Base/Dijkstra.hxx
new file mode 100644 (file)
index 0000000..4920ab5
--- /dev/null
@@ -0,0 +1,95 @@
+#ifndef __FPA__BASE__DIJKSTRA__HXX__
+#define __FPA__BASE__DIJKSTRA__HXX__
+
+#include <algorithm>
+
+// -------------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+fpa::Base::Dijkstra< V, C, VV, VC, B >::
+Dijkstra( )
+  : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+fpa::Base::Dijkstra< V, C, VV, VC, B >::
+~Dijkstra( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+void fpa::Base::Dijkstra< V, C, VV, VC, B >::
+_InitializeQueue( )
+{
+  for(
+    typename _TNodes::const_iterator vIt = this->m_Seeds.begin( );
+    vIt != this->m_Seeds.end( );
+    vIt++
+    )
+    this->_QueuePush( *vIt );
+}
+
+// -------------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+bool fpa::Base::Dijkstra< V, C, VV, VC, B >::
+_IsQueueEmpty( ) const
+{
+  return( this->m_Queue.empty( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+void fpa::Base::Dijkstra< V, C, VV, VC, B >::
+_QueuePush( const _TNode& n )
+{
+  this->m_Queue.push_back( n );
+  std::push_heap( this->m_Queue.begin( ), this->m_Queue.end( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+typename fpa::Base::Dijkstra< V, C, VV, VC, B >::
+_TNode fpa::Base::Dijkstra< V, C, VV, VC, B >::
+_QueuePop( )
+{
+  _TNode n;
+  if( !( this->m_Queue.empty( ) ) )
+  {
+    // n = *( this->m_Queue.begin( ) );
+    n = this->m_Queue.front( );
+    std::pop_heap( this->m_Queue.begin( ), this->m_Queue.end( ) );
+    this->m_Queue.pop_back( );
+
+  } // fi
+  return( n );
+}
+
+// -------------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+void fpa::Base::Dijkstra< V, C, VV, VC, B >::
+_QueueClear( )
+{
+  this->m_Queue.clear( );
+}
+
+// -------------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+bool fpa::Base::Dijkstra< V, C, VV, VC, B >::
+_UpdateNeigh( _TNode& nn, const _TNode& n )
+{
+  TCost nc = this->_Cost( nn.Vertex, n.Vertex );
+  if( C( 0 ) <= nc )
+  {
+    nn.Cost = n.Cost + nc;
+    nn.Result = nn.Cost;
+    return( true );
+  }
+  else
+    return( false );
+}
+
+#endif // __FPA__BASE__DIJKSTRA__HXX__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Base/Events.h b/lib/fpa/Base/Events.h
new file mode 100644 (file)
index 0000000..b7b9f4e
--- /dev/null
@@ -0,0 +1,136 @@
+#ifndef __FPA__BASE__EVENTS__H__
+#define __FPA__BASE__EVENTS__H__
+
+#include <vector>
+#include <itkProcessObject.h>
+
+namespace fpa
+{
+  namespace Base
+  {
+    /**
+     * Evolution event. An event is generated when a vertex changes its
+     * state.
+     */
+    template< class N >
+    class BaseEvent
+      : public itk::AnyEvent
+    {
+    public:
+      BaseEvent( )
+        : itk::AnyEvent( )
+        { }
+      BaseEvent( const N& n )
+        : itk::AnyEvent( ),
+          Node( n )
+        { }
+      virtual ~BaseEvent( )
+        { }
+
+      const char* GetEventName( ) const
+        { return( "fpa::Base::BaseEvent" ); }
+      bool CheckEvent( const itk::EventObject* e ) const
+        { return( dynamic_cast< const BaseEvent< N >* >( e ) != NULL ); }
+      itk::EventObject* MakeObject( ) const
+        { return( new BaseEvent< N >( ) ); }
+
+    public:
+      N Node;
+    };
+
+    /**
+     */
+    template< class N >
+    class FrontEvent
+      : public BaseEvent< N >
+    {
+    public:
+      FrontEvent( )
+        : BaseEvent< N >( )
+        { }
+      FrontEvent( const N& n )
+        : BaseEvent< N >( n )
+        { }
+      virtual ~FrontEvent( )
+        { }
+      const char* GetEventName( ) const
+        { return( "fpa::Base::FrontEvent" ); }
+      bool CheckEvent( const itk::EventObject* e ) const
+        { return( dynamic_cast< const FrontEvent< N >* >( e ) != NULL ); }
+      itk::EventObject* MakeObject( ) const
+        { return( new FrontEvent< N >( ) ); }
+    };
+
+    /**
+     */
+    template< class N >
+    class MarkEvent
+      : public BaseEvent< N >
+    {
+    public:
+      MarkEvent( )
+        : BaseEvent< N >( )
+        { }
+      MarkEvent( const N& n )
+        : BaseEvent< N >( n )
+        { }
+      virtual ~MarkEvent( )
+        { }
+      const char* GetEventName( ) const
+        { return( "fpa::Base::MarkEvent" ); }
+      bool CheckEvent( const itk::EventObject* e ) const
+        { return( dynamic_cast< const MarkEvent< N >* >( e ) != NULL ); }
+      itk::EventObject* MakeObject( ) const
+        { return( new MarkEvent< N >( ) ); }
+    };
+
+    /**
+     */
+    template< class N >
+    class CollisionEvent
+      : public BaseEvent< N >
+    {
+    public:
+      CollisionEvent( )
+        : BaseEvent< N >( )
+        { }
+      CollisionEvent( const N& n )
+        : BaseEvent< N >( n )
+        { }
+      virtual ~CollisionEvent( )
+        { }
+      const char* GetEventName( ) const
+        { return( "fpa::Base::CollisionEvent" ); }
+      bool CheckEvent( const itk::EventObject* e ) const
+        { return( dynamic_cast< const CollisionEvent< N >* >( e ) != NULL ); }
+      itk::EventObject* MakeObject( ) const
+        { return( new CollisionEvent< N >( ) ); }
+    };
+
+    /**
+     */
+    template< class N >
+    class EndEvent
+      : public BaseEvent< N >
+    {
+    public:
+      EndEvent( )
+        : BaseEvent< N >( )
+        { }
+      virtual ~EndEvent( )
+        { }
+      const char* GetEventName( ) const
+        { return( "fpa::Base::EndEvent" ); }
+      bool CheckEvent( const itk::EventObject* e ) const
+        { return( dynamic_cast< const EndEvent< N >* >( e ) != NULL ); }
+      itk::EventObject* MakeObject( ) const
+        { return( new EndEvent< N >( ) ); }
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#endif // __FPA__BASE__EVENTS__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Base/FastMarching.h b/lib/fpa/Base/FastMarching.h
new file mode 100644 (file)
index 0000000..119d0f3
--- /dev/null
@@ -0,0 +1,76 @@
+#ifndef __FPA__BASE__FASTMARCHING__H__
+#define __FPA__BASE__FASTMARCHING__H__
+
+#include <fpa/Base/Dijkstra.h>
+
+namespace fpa
+{
+  namespace Base
+  {
+    /**
+     */
+    template< class V, class C, class VV, class VC, class B >
+    class FastMarching
+      : public Dijkstra< V, C, VV, VC, B >
+    {
+    public:
+      // Templated types
+      typedef V  TVertex;
+      typedef C  TCost;
+      typedef VV TVertexValue;
+      typedef B  TBaseFilter;
+
+      // Standard class typdedefs
+      typedef FastMarching                    Self;
+      typedef Dijkstra< V, C, VV, VC, B >     Superclass;
+      typedef itk::SmartPointer< Self >       Pointer;
+      typedef itk::SmartPointer< const Self > ConstPointer;
+
+      typedef typename Superclass::TResult TResult;
+      typedef typename Superclass::TTraits TTraits;
+
+    protected:
+      typedef typename Superclass::_TFrontId _TFrontId;
+      typedef typename Superclass::_TNode    _TNode;
+      typedef typename Superclass::_TNodes   _TNodes;
+      struct _TNodeEikonalSort
+      {
+        bool operator()( const _TNode& a, const _TNode& b ) const
+          { return( a.Result < b.Result ); }
+      };
+
+    public:
+      itkTypeMacro( FastMarching, Dijkstra );
+
+    protected:
+      FastMarching( );
+      virtual ~FastMarching( );
+
+      virtual bool _UpdateNeigh( _TNode& nn, const _TNode& n );
+
+      /**
+       * Eikonal solution method
+       *
+       * This method returns the computed time and a flag indicating if the
+       * determinant was valid.
+       */
+      bool _SolveEikonal( double& re, const _TNode& n, const double& F );
+
+    private:
+      // Purposely not implemented
+      FastMarching( const Self& );
+      void operator=( const Self& );
+
+    protected:
+      static const double LargestValue;
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#include <fpa/Base/FastMarching.hxx>
+
+#endif // __FPA__BASE__FASTMARCHING__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Base/FastMarching.hxx b/lib/fpa/Base/FastMarching.hxx
new file mode 100644 (file)
index 0000000..ef1f634
--- /dev/null
@@ -0,0 +1,121 @@
+#ifndef __FPA__BASE__FASTMARCHING__HXX__
+#define __FPA__BASE__FASTMARCHING__HXX__
+
+#include <cmath>
+#include <limits>
+
+#include <itkExceptionObject.h>
+
+// ---------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+const double fpa::Base::FastMarching< V, C, VV, VC, B >::
+LargestValue = std::numeric_limits< double >::max( );
+
+// -------------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+fpa::Base::FastMarching< V, C, VV, VC, B >::
+FastMarching( )
+  : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+fpa::Base::FastMarching< V, C, VV, VC, B >::
+~FastMarching( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+bool fpa::Base::FastMarching< V, C, VV, VC, B >::
+_UpdateNeigh( _TNode& nn, const _TNode& n )
+{
+  bool r = true;
+  if( !( this->_IsMarked( nn.Vertex ) ) )
+  {
+    double re;
+    double F = double( this->_Value( nn.Vertex ) );
+    r = this->_SolveEikonal( re, nn, F );
+    if( r )
+    {
+      nn.Result = TResult( re );
+      nn.Cost = TCost( re );
+
+    } // fi
+
+  } // fi
+  return( r );
+}
+
+// -------------------------------------------------------------------------
+template< class V, class C, class VV, class VC, class B >
+bool fpa::Base::FastMarching< V, C, VV, VC, B >::
+_SolveEikonal( double& re, const _TNode& n, const double& F )
+{
+  _TNodes neighs;
+  this->_Neighs( n, neighs );
+
+  // Compute interest nodes (i.e. the smallest valued in each dimension)
+  _TNodes minNodes;
+  unsigned int i = 0;
+  while( i < neighs.size( ) )
+  {
+    V nm = neighs[ i++ ].Vertex;
+    V np = neighs[ i++ ].Vertex;
+    bool mm = this->_IsMarked( nm );
+    bool mp = this->_IsMarked( np );
+    TResult rm = this->_Result( nm );
+    TResult rp = this->_Result( np );
+
+    _TNode Ti;
+    if( mm && mp )
+      minNodes.push_back(
+        ( rm < rp )?
+        _TNode( nm, rm, n.FrontId ):
+        _TNode( np, rp, n.FrontId )
+        );
+    else if( !mm && mp )
+      minNodes.push_back( _TNode( np, rp, n.FrontId ) );
+    else if( mm && !mp )
+      minNodes.push_back( _TNode( nm, rm, n.FrontId ) );
+
+  } // elihw
+  std::sort( minNodes.begin( ), minNodes.end( ), _TNodeEikonalSort( ) );
+
+  // Compute solution
+  double al = double( 0 );
+  double be = double( 0 );
+  double ga = double( -1 ) / ( F * F );
+  re = Self::LargestValue;
+
+  typename _TNodes::const_iterator nIt = minNodes.begin( );
+  bool valid = true;
+  while( nIt != minNodes.end( ) && valid )
+  {
+    if( !( re < double( nIt->Result ) ) )
+    {
+      double d = this->_Norm( n.Vertex, nIt->Vertex );
+      double s = double( 1 ) / ( d * d );
+      double v = double( nIt->Result );
+
+      al += s;
+      be += v * s;
+      ga += v * v * s;
+
+      double de = ( be * be ) - ( al * ga );
+      valid = !( de < double( 0 ) );
+      if( valid )
+        re = ( be + std::sqrt( de ) ) / al;
+      nIt++;
+    }
+    else
+      nIt = minNodes.end( );
+
+  } // elihw
+  return( valid );
+}
+
+#endif // __FPA__BASE__FASTMARCHING__HXX__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Base/Functors/InvertCostFunction.h b/lib/fpa/Base/Functors/InvertCostFunction.h
new file mode 100644 (file)
index 0000000..faf4fb2
--- /dev/null
@@ -0,0 +1,56 @@
+#ifndef __FPA__BASE__FUNCTORS__INVERTCOSTFUNCTION__H__
+#define __FPA__BASE__FUNCTORS__INVERTCOSTFUNCTION__H__
+
+#include <itkFunctionBase.h>
+
+namespace fpa
+{
+  namespace Base
+  {
+    namespace Functors
+    {
+      /**
+       */
+      template< class C >
+      class InvertCostFunction
+        : public itk::FunctionBase< C, C >
+      {
+      public:
+        /// Type-related and pointers
+        typedef InvertCostFunction              Self;
+        typedef itk::FunctionBase< C, C >       Superclass;
+        typedef itk::SmartPointer< Self >       Pointer;
+        typedef itk::SmartPointer< const Self > ConstPointer;
+
+      public:
+        itkNewMacro( Self );
+        itkTypeMacro( InvertCostFunction, itkFunctionBase );
+
+      public:
+        virtual C Evaluate( const C& input ) const
+          {
+            return( C( 1 ) / ( C( 1 ) + C( input ) ) );
+          }
+
+      protected:
+        InvertCostFunction( )
+          : Superclass( )
+          { }
+        virtual ~InvertCostFunction( )
+          { }
+
+      private:
+        // Purposely not implemented
+        InvertCostFunction( const Self& );
+        void operator=( const Self& );
+      };
+
+    } // ecapseman
+
+  } // ecapseman
+
+} // ecapseman
+
+#endif // __FPA__BASE__FUNCTORS__INVERTCOSTFUNCTION__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Base/RegionGrow.h b/lib/fpa/Base/RegionGrow.h
new file mode 100644 (file)
index 0000000..e52fb62
--- /dev/null
@@ -0,0 +1,117 @@
+#ifndef __FPA__BASE__REGIONGROW__H__
+#define __FPA__BASE__REGIONGROW__H__
+
+#include <queue>
+#include <vector>
+#include <fpa/Base/Algorithm.h>
+
+namespace fpa
+{
+  namespace Base
+  {
+    /**
+     */
+    template< class V, class R, class VV, class VC >
+    class RegionGrowTraits
+    {
+    public:
+      typedef R  TResult;
+      typedef V  TVertex;
+      typedef VV TVertexValue;
+      typedef VC TVertexCmp;
+
+      typedef bool TCost;
+      typedef long TFrontId;
+
+      class TNode
+      {
+      public:
+        TNode( )
+          { }
+        TNode( const TVertex& v, const TFrontId& f )
+          : Vertex( v ),
+            Parent( v ),
+            FrontId( f )
+          { }
+        TNode( const TVertex& v, const TResult& r, const TFrontId& f )
+          : Vertex( v ),
+            Parent( v ),
+            Result( r ),
+            FrontId( f )
+          { }
+        virtual ~TNode( )
+          { }
+
+        TVertex Vertex;
+        TVertex Parent;
+        TResult Result;
+        TFrontId FrontId;
+      };
+
+      typedef std::vector< TNode > TNodes;
+    };
+
+    /**
+     * Region grow is a front propagation with no costs.
+     */
+    template< class V, class R, class VV, class VC, class B >
+    class RegionGrow
+      : public Algorithm< RegionGrowTraits< V, R, VV, VC >, B >
+    {
+    public:
+      typedef V  TVertex;
+      typedef R  TResult;
+      typedef VV TVertexValue;
+      typedef B  TBaseFilter;
+
+      /// Standard class typdedefs
+      typedef RegionGrowTraits< V, R, VV, VC > TTraits;
+      typedef RegionGrow                       Self;
+      typedef Algorithm< TTraits, B >          Superclass;
+      typedef itk::SmartPointer< Self >        Pointer;
+      typedef itk::SmartPointer< const Self >  ConstPointer;
+
+      typedef typename TTraits::TCost TCost;
+
+    protected:
+      typedef typename TTraits::TFrontId  _TFrontId;
+      typedef typename TTraits::TNode     _TNode;
+      typedef typename TTraits::TNodes    _TNodes;
+
+    private:
+      typedef std::queue< _TNode > _TQueue;
+
+    public:
+      itkTypeMacro( RegionGrow, Base );
+
+    protected:
+      RegionGrow( );
+      virtual ~RegionGrow( );
+
+      virtual   bool _UpdateResult    ( _TNode& n );
+      virtual   void _InitializeQueue ( );
+      virtual   bool _IsQueueEmpty    ( ) const;
+      virtual   void _QueuePush       ( const _TNode& n );
+      virtual _TNode _QueuePop        ( );
+      virtual   void _QueueClear      ( );
+      virtual   bool _UpdateNeigh     ( _TNode& nn, const _TNode& n );
+
+      virtual bool _CheckMembership( const _TNode& n ) const = 0;
+
+    private:
+      RegionGrow( const Self& );     // Not impl.
+      void operator=( const Self& ); // Not impl.
+
+    private:
+      _TQueue m_Queue;
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#include <fpa/Base/RegionGrow.hxx>
+
+#endif // __FPA__BASE__REGIONGROW__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Base/RegionGrow.hxx b/lib/fpa/Base/RegionGrow.hxx
new file mode 100644 (file)
index 0000000..bdc1ed3
--- /dev/null
@@ -0,0 +1,87 @@
+#ifndef __FPA__BASE__REGIONGROWING__HXX__
+#define __FPA__BASE__REGIONGROWING__HXX__
+
+// -------------------------------------------------------------------------
+template< class V, class R, class VV, class VC, class B >
+fpa::Base::RegionGrow< V, R, VV, VC, B >::
+RegionGrow( )
+  : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class V, class R, class VV, class VC, class B >
+fpa::Base::RegionGrow< V, R, VV, VC, B >::
+~RegionGrow( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class V, class R, class VV, class VC, class B >
+bool fpa::Base::RegionGrow< V, R, VV, VC, B >::
+_UpdateResult( _TNode& n )
+{
+  n.Result = R( this->_CheckMembership( n ) );
+  return( n.Result );
+}
+
+// -------------------------------------------------------------------------
+template< class V, class R, class VV, class VC, class B >
+void fpa::Base::RegionGrow< V, R, VV, VC, B >::
+_InitializeQueue( )
+{
+  for(
+    typename _TNodes::const_iterator vIt = this->m_Seeds.begin( );
+    vIt != this->m_Seeds.end( );
+    vIt++
+    )
+    this->_QueuePush( *vIt );
+}
+
+// -------------------------------------------------------------------------
+template< class V, class R, class VV, class VC, class B >
+bool fpa::Base::RegionGrow< V, R, VV, VC, B >::
+_IsQueueEmpty( ) const
+{
+  return( this->m_Queue.empty( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class V, class R, class VV, class VC, class B >
+void fpa::Base::RegionGrow< V, R, VV, VC, B >::
+_QueuePush( const _TNode& n )
+{
+  this->m_Queue.push( n );
+}
+
+// -------------------------------------------------------------------------
+template< class V, class R, class VV, class VC, class B >
+typename fpa::Base::RegionGrow< V, R, VV, VC, B >::
+_TNode fpa::Base::RegionGrow< V, R, VV, VC, B >::
+_QueuePop( )
+{
+  _TNode n = this->m_Queue.front( );
+  this->m_Queue.pop( );
+  return( n );
+}
+
+// -------------------------------------------------------------------------
+template< class V, class R, class VV, class VC, class B >
+void fpa::Base::RegionGrow< V, R, VV, VC, B >::
+_QueueClear( )
+{
+  while( this->m_Queue.size( ) > 0 )
+    this->m_Queue.pop( );
+}
+
+// -------------------------------------------------------------------------
+template< class V, class R, class VV, class VC, class B >
+bool fpa::Base::RegionGrow< V, R, VV, VC, B >::
+_UpdateNeigh( _TNode& nn, const _TNode& n )
+{
+  return( true );
+}
+
+#endif // __FPA__BASE__REGIONGROWING__HXX__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Base/TreeExtractor.h b/lib/fpa/Base/TreeExtractor.h
new file mode 100644 (file)
index 0000000..a0222a6
--- /dev/null
@@ -0,0 +1,84 @@
+#ifndef __FPA__BASE__TREEEXTRACTOR__H__
+#define __FPA__BASE__TREEEXTRACTOR__H__
+
+#include <map>
+#include <vector>
+
+namespace fpa
+{
+  namespace Base
+  {
+    /**
+     */
+    template< class A >
+    class TreeExtractor
+      : public A
+    {
+    public:
+      /// Standard class typdedefs
+      typedef TreeExtractor                   Self;
+      typedef A                               Superclass;
+      typedef itk::SmartPointer< Self >       Pointer;
+      typedef itk::SmartPointer< const Self > ConstPointer;
+
+      typedef A TAlgorithm;
+
+      typedef typename A::TCost                        TCost;
+      typedef typename A::TVertex                      TVertex;
+      typedef std::vector< TVertex >                   TVertices;
+      typedef typename A::TTraits::TVertexCmp          TVertexCmp;
+      typedef std::map< TVertex, TVertex, TVertexCmp > TTree;
+
+    protected:
+      typedef typename A::_TFrontId           _TFrontId;
+      typedef typename A::_TNode              _TNode;
+      typedef typename A::_TCollisionSitesRow _TCollisionSitesRow;
+      typedef typename A::_TCollisionSites    _TCollisionSites;
+      typedef std::vector< unsigned long >    _TRow;
+      typedef std::vector< _TRow >            _TMatrix;
+
+    public:
+      itkNewMacro( Self );
+      itkTypeMacro( TreeExtractor, TAlgorithm );
+
+      itkGetConstMacro( Root, TVertex );
+      itkGetConstMacro( Tree, TTree );
+      itkSetMacro( Root, TVertex );
+
+    public:
+      unsigned int AddLeaf( const TVertex& leaf );
+      unsigned int GetNumberOfLeafs( ) const;
+      const TVertex& GetLeaf( const unsigned int& id ) const;
+
+    protected:
+      TreeExtractor( );
+      virtual ~TreeExtractor( );
+
+      virtual void _BeforeLoop( );
+      virtual void _AfterLoop( );
+
+      TVertices _Path( const TVertex& s, const TVertex& e );
+
+    private:
+      // Purporsely not implemented
+      TreeExtractor( const Self& );
+      void operator=( const Self& );
+
+    protected:
+      TVertex   m_Root;
+      TVertices m_Leafs;
+      TTree     m_Tree;
+
+      _TMatrix m_FrontPaths;
+      static const unsigned long INF_VALUE;
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#include <fpa/Base/TreeExtractor.hxx>
+
+#endif // __FPA__BASE__TREEEXTRACTOR__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Base/TreeExtractor.hxx b/lib/fpa/Base/TreeExtractor.hxx
new file mode 100644 (file)
index 0000000..cfeb17f
--- /dev/null
@@ -0,0 +1,250 @@
+#ifndef __FPA__BASE__TREEEXTRACTOR__HXX__
+#define __FPA__BASE__TREEEXTRACTOR__HXX__
+
+#include <limits>
+
+// -------------------------------------------------------------------------
+template< class A >
+const unsigned long fpa::Base::TreeExtractor< A >::INF_VALUE =
+  std::numeric_limits< unsigned long >::max( ) >> 1;
+
+// -------------------------------------------------------------------------
+template< class A >
+unsigned int fpa::Base::TreeExtractor< A >::
+AddLeaf( const TVertex& leaf )
+{
+  this->m_Leafs.push_back( leaf );
+  this->Modified( );
+  return( this->m_Leafs.size( ) - 1 );
+}
+
+// -------------------------------------------------------------------------
+template< class A >
+unsigned int fpa::Base::TreeExtractor< A >::
+GetNumberOfLeafs( ) const
+{
+  return( this->m_Leafs.size( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class A >
+const typename fpa::Base::TreeExtractor< A >::
+TVertex& fpa::Base::TreeExtractor< A >::
+GetLeaf( const unsigned int& id ) const
+{
+  static const TVertex zero;
+  if( id < this->m_Leafs.size( ) )
+    return( this->m_Leafs[ id ] );
+  else
+    return( zero );
+}
+
+// -------------------------------------------------------------------------
+template< class A >
+fpa::Base::TreeExtractor< A >::
+TreeExtractor( )
+  : Superclass( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class A >
+fpa::Base::TreeExtractor< A >::
+~TreeExtractor( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class A >
+void fpa::Base::TreeExtractor< A >::
+_BeforeLoop( )
+{
+  this->Superclass::_BeforeLoop( );
+  this->m_Tree.clear( );
+}
+
+// -------------------------------------------------------------------------
+template< class A >
+void fpa::Base::TreeExtractor< A >::
+_AfterLoop( )
+{
+  this->Superclass::_AfterLoop( );
+  
+  // Prepare fronts graph using Floyd-Warshall
+  unsigned long nSeeds = this->GetNumberOfSeeds( );
+  _TMatrix dist( nSeeds, _TRow( nSeeds, Self::INF_VALUE ) );
+  this->m_FrontPaths = dist;
+
+  for( unsigned long i = 0; i < nSeeds; ++i )
+  {
+    for( unsigned long j = 0; j < nSeeds; ++j )
+    {
+      if( this->m_CollisionSites[ i ][ j ].second )
+      {
+        dist[ i ][ j ] = 1;
+        dist[ j ][ i ] = 1;
+        this->m_FrontPaths[ i ][ j ] = j;
+        this->m_FrontPaths[ j ][ i ] = i;
+
+      } // fi
+
+    } // rof
+    dist[ i ][ i ] = 0;
+    this->m_FrontPaths[ i ][ i ] = i;
+
+  } // rof
+  for( unsigned long k = 0; k < nSeeds; ++k )
+  {
+    for( unsigned long i = 0; i < nSeeds; ++i )
+    {
+      for( unsigned long j = 0; j < nSeeds; ++j )
+      {
+        if( ( dist[ i ][ k ] + dist[ k ][ j ] ) < dist[ i ][ j ] )
+        {
+          dist[ i ][ j ] = dist[ i ][ k ] + dist[ k ][ j ];
+          this->m_FrontPaths[ i ][ j ] = this->m_FrontPaths[ i ][ k ];
+
+        } // fi
+
+      } // rof
+
+    } // rof
+
+  } // rof
+
+  // Build tree
+  this->m_Tree.clear( );
+  typename TVertices::const_iterator vIt = this->m_Leafs.begin( );
+  for( ; vIt != this->m_Leafs.end( ); ++vIt )
+  {
+    TVertices path = this->_Path( this->m_Root, *vIt );
+    typename TVertices::const_iterator pIt = path.begin( );
+    TVertex parent = *pIt;
+    for( ; pIt != path.end( ); ++pIt )
+    {
+      this->m_Tree[ *pIt ] = parent;
+      parent = *pIt;
+
+    } // rof
+
+  } // rof
+}
+
+// -------------------------------------------------------------------------
+template< class A >
+typename fpa::Base::TreeExtractor< A >::
+TVertices fpa::Base::TreeExtractor< A >::
+_Path( const TVertex& s, const TVertex& e )
+{
+  TVertices pS, pE;
+
+  // Check if both vertices belong to the same front or if collisions
+  // should be analyzed.
+  _TFrontId fs = this->_FrontId( s );
+  _TFrontId fe = this->_FrontId( e );
+  if( fs == fe )
+  {
+    // Path from start vertex until seed
+    TVertex vIt = s;
+    TVertex parent = this->_Parent( vIt );
+    while( parent != vIt )
+    {
+      pS.push_back( vIt );
+      vIt = parent;
+      parent = this->_Parent( vIt );
+
+    } // elihw
+    pS.push_back( vIt );
+
+    // Path from end vertex until seed
+    vIt = e;
+    parent = this->_Parent( vIt );
+    while( parent != vIt )
+    {
+      pE.push_back( vIt );
+      vIt = parent;
+      parent = this->_Parent( vIt );
+
+    } // elihw
+    pE.push_back( vIt );
+
+    // Erase duplicate vertices
+    if( pS.size( ) > 0 && pE.size( ) > 0 )
+    {
+      bool cont = ( pS.back( ) == pE.back( ) );
+      bool enter = false;
+      TVertex lastVertex = pS.back( );
+      while( cont )
+      {
+        enter = true;
+
+        lastVertex = pS.back( );
+        pS.pop_back( );
+        pE.pop_back( );
+
+        cont = ( pS.size( ) > 0 && pE.size( ) > 0 );
+        if( cont )
+          cont = ( pS.back( ) == pE.back( ) );
+
+      } // elihw
+      if( enter )
+        pS.push_back( lastVertex );
+
+    } // fi
+
+    // Contatenate both paths
+    pS.insert( pS.end( ), pE.rbegin( ), pE.rend( ) );
+  }
+  else
+  {
+    // Use this->m_FrontPaths from Floyd-Warshall
+    if( this->m_FrontPaths[ fs ][ fe ] != Self::INF_VALUE )
+    {
+      // Compute front path
+      std::vector< unsigned long > fpath;
+      fpath.push_back( fs );
+      while( fs != fe )
+      {
+        fs = this->m_FrontPaths[ fs ][ fe ];
+        fpath.push_back( fs );
+
+      } // elihw
+
+      // Continue only if both fronts are connected
+      unsigned int N = fpath.size( );
+      if( 0 < N )
+      {
+        // First path: from start vertex to first collision
+        pS = this->_Path(
+          s, this->m_CollisionSites[ fpath[ 0 ] ][ fpath[ 1 ] ].first
+          );
+
+        // Intermediary paths
+        for( unsigned int i = 1; i < N - 1; ++i )
+        {
+          pE = this->_Path(
+            this->m_CollisionSites[ fpath[ i ] ][ fpath[ i - 1 ] ].first,
+            this->m_CollisionSites[ fpath[ i ] ][ fpath[ i + 1 ] ].first
+            );
+          pS.insert( pS.end( ), pE.begin( ), pE.end( ) );
+
+        } // rof
+
+        // Final path: from last collision to end point
+        pE = this->_Path(
+          this->m_CollisionSites[ fpath[ N - 1 ] ][ fpath[ N - 2 ] ].first,
+          e
+          );
+        pS.insert( pS.end( ), pE.begin( ), pE.end( ) );
+
+      } // fi
+
+    } // fi
+
+  } // fi
+  return( pS );
+}
+
+#endif // __FPA__BASE__TREEEXTRACTOR__HXX__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Common.cxx.in b/lib/fpa/Common.cxx.in
new file mode 100644 (file)
index 0000000..2b1a193
--- /dev/null
@@ -0,0 +1,9 @@
+#include <fpa/Common.h>
+
+// -------------------------------------------------------------------------
+std::string fpa::GetVersion( )
+{
+  return( "@VERSION_STRING@" );
+}
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Common.h b/lib/fpa/Common.h
new file mode 100644 (file)
index 0000000..c5bf168
--- /dev/null
@@ -0,0 +1,15 @@
+#ifndef __FPA__COMMON__H__
+#define __FPA__COMMON__H__
+
+#include <string>
+#include <fpa/FrontAlgorithms_Export.h>
+
+namespace fpa
+{
+  std::string FrontAlgorithms_EXPORT GetVersion( );
+
+} // ecapseman
+
+#endif // __FPA__COMMON__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/Algorithm.h b/lib/fpa/Image/Algorithm.h
new file mode 100644 (file)
index 0000000..289eb6c
--- /dev/null
@@ -0,0 +1,105 @@
+#ifndef __FPA__IMAGE__ALGORITHM__H__
+#define __FPA__IMAGE__ALGORITHM__H__
+
+#include <itkImage.h>
+#include <itkFunctionBase.h>
+
+namespace fpa
+{
+  namespace Image
+  {
+    /**
+     * A generic front propagation algorithm were vertices are image pixels.
+     *
+     * @param I    Input image type
+     * @param A    Base algorithm (RegionGrow, Dijkstra or FastMarching)
+     */
+    template< class I, class A >
+    class Algorithm
+      : public A
+    {
+    public:
+      /// Standard class typdedefs
+      typedef Algorithm                       Self;
+      typedef A                               Superclass;
+      typedef itk::SmartPointer< Self >       Pointer;
+      typedef itk::SmartPointer< const Self > ConstPointer;
+
+      /// Template input values
+      typedef I TInputImage;
+      typedef A TBaseAlgorithm;
+
+      typedef typename A::TTraits            TTraits;
+      typedef typename TTraits::TCost        TCost;
+      typedef typename TTraits::TResult      TResult;
+      typedef typename TTraits::TVertex      TVertex;
+      typedef typename TTraits::TVertexValue TVertexValue;
+
+      typedef itk::Image< TResult, I::ImageDimension > TOutputImage;
+      typedef itk::FunctionBase< TCost, TCost > TCostConversionFunction;
+
+    protected:
+      typedef typename TTraits::TFrontId _TFrontId;
+      typedef typename TTraits::TNode    _TNode;
+      typedef typename TTraits::TNodes   _TNodes;
+
+    private:
+      typedef itk::Image< bool, I::ImageDimension >      _TMarks;
+      typedef itk::Image< _TFrontId, I::ImageDimension > _TFrontsIds;
+      typedef itk::Image< TVertex, I::ImageDimension >   _TParents;
+
+    public:
+      itkTypeMacro( Algorithm, TAlgorithm );
+
+      /// Set/Get
+      itkGetConstMacro( NeighborhoodOrder, unsigned int );
+      itkGetConstObjectMacro( CostConversion, TCostConversionFunction );
+      itkGetObjectMacro( CostConversion, TCostConversionFunction );
+
+      itkSetMacro( NeighborhoodOrder, unsigned int );
+      itkSetObjectMacro( CostConversion, TCostConversionFunction );
+
+    protected:
+      Algorithm( );
+      virtual ~Algorithm( );
+
+      /// Base interface
+      virtual bool _UpdateResult( _TNode& n );
+
+      /// Pure virtual interface: vertices
+      virtual unsigned long _NumberOfVertices ( ) const;
+      virtual  TVertexValue _Value            ( const TVertex& v ) const;
+      virtual       TResult _Result           ( const TVertex& v ) const;
+
+      /// Pure virtual interface: edges
+      virtual double _Norm ( const TVertex& a, const TVertex& b ) const;
+      virtual   bool _Edge ( const TVertex& a, const TVertex& b ) const;
+      virtual  TCost _Cost ( const TVertex& a, const TVertex& b ) const;
+
+      /// Pure virtual interface: neighborhood
+      virtual void _Neighs      ( const _TNode& n, _TNodes& N ) const;
+      virtual void _NeighsInDim ( const _TNode& n,
+                                  const unsigned int& d,
+                                  _TNodes& N );
+
+      /// Pure virtual interface: results
+      virtual void _InitializeResults ( );
+
+    private:
+      Algorithm( const Self& );      // Not impl.
+      void operator=( const Self& ); // Not impl.
+
+    protected:
+      unsigned int                              m_NeighborhoodOrder;
+      typename TCostConversionFunction::Pointer m_CostConversion;
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#include <fpa/Image/Algorithm.hxx>
+
+#endif // __FPA__IMAGE__ALGORITHM__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/Algorithm.hxx b/lib/fpa/Image/Algorithm.hxx
new file mode 100644 (file)
index 0000000..7800a54
--- /dev/null
@@ -0,0 +1,180 @@
+#ifndef __FPA__IMAGE__ALGORITHM__HXX__
+#define __FPA__IMAGE__ALGORITHM__HXX__
+
+#include <cmath>
+#include <limits>
+#include <itkConstNeighborhoodIterator.h>
+
+// -------------------------------------------------------------------------
+template< class I, class A >
+fpa::Image::Algorithm< I, A >::
+Algorithm( )
+  : Superclass( ),
+    m_NeighborhoodOrder( 1 )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class I, class A >
+fpa::Image::Algorithm< I, A >::
+~Algorithm( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class I, class A >
+bool fpa::Image::Algorithm< I, A >::
+_UpdateResult( _TNode& n )
+{
+  bool ret = this->Superclass::_UpdateResult( n );
+  this->GetOutput( )->SetPixel( n.Vertex, n.Result );
+  return( ret );
+}
+
+// -------------------------------------------------------------------------
+template< class I, class A >
+unsigned long fpa::Image::Algorithm< I, A >::
+_NumberOfVertices( ) const
+{
+  return(
+    this->GetInput( )->GetLargestPossibleRegion( ).GetNumberOfPixels( )
+    );
+}
+
+// -------------------------------------------------------------------------
+template< class I, class A >
+typename fpa::Image::Algorithm< I, A >::
+TVertexValue fpa::Image::Algorithm< I, A >::
+_Value( const TVertex& v ) const
+{
+  return( this->GetInput( )->GetPixel( v ) );
+}
+
+// -------------------------------------------------------------------------
+template< class I, class A >
+typename fpa::Image::Algorithm< I, A >::
+TResult fpa::Image::Algorithm< I, A >::
+_Result( const TVertex& v ) const
+{
+  return( this->GetOutput( )->GetPixel( v ) );
+}
+
+// -------------------------------------------------------------------------
+template< class I, class A >
+double fpa::Image::Algorithm< I, A >::
+_Norm( const TVertex& a, const TVertex& b ) const
+{
+  typename I::PointType pa, pb;
+  this->GetInput( )->TransformIndexToPhysicalPoint( a, pa );
+  this->GetInput( )->TransformIndexToPhysicalPoint( b, pb );
+  return( double( pa.EuclideanDistanceTo( pb ) ) );
+}
+
+// -------------------------------------------------------------------------
+template< class I, class A >
+bool fpa::Image::Algorithm< I, A >::
+_Edge( const TVertex& a, const TVertex& b ) const
+{
+  unsigned long dist = 0;
+  for( unsigned int d = 0; d < I::ImageDimension; d++ )
+    dist += std::abs( long( a[ d ] ) - long( b[ d ] ) );
+  if( this->m_NeighborhoodOrder == 1 )
+    return( dist == 0 || dist == 1 );
+  else
+    return( dist == 0 || dist == 1 || dist == 2 );
+}
+
+// -------------------------------------------------------------------------
+template< class I, class A >
+typename fpa::Image::Algorithm< I, A >::
+TCost fpa::Image::Algorithm< I, A >::
+_Cost( const TVertex& a, const TVertex& b ) const
+{
+  static const TCost INF_COST = std::numeric_limits< TCost >::max( );
+  if( this->_Edge( a, b ) )
+  {
+    TCost c = TCost( this->GetInput( )->GetPixel( b ) );
+    if( this->m_CostConversion.IsNotNull( ) )
+      return( this->m_CostConversion->Evaluate( c ) );
+    else
+      return( c );
+  }
+  else
+    return( INF_COST );
+}
+
+// -------------------------------------------------------------------------
+template< class I, class A >
+void fpa::Image::Algorithm< I, A >::
+_Neighs( const _TNode& n, _TNodes& N ) const
+{
+  typename I::RegionType reg = this->GetInput( )->GetRequestedRegion( );
+
+  N.clear( );
+  if( this->m_NeighborhoodOrder == 1 )
+  {
+    for( unsigned int d = 0; d < I::ImageDimension; d++ )
+    {
+      for( int i = -1; i <= 1; i += 2 )
+      {
+        TVertex v = n.Vertex;
+        v[ d ] += i;
+        if( reg.IsInside( v ) )
+          N.push_back( _TNode( v, n.FrontId ) );
+        else
+          N.push_back( n );
+
+      } // rof
+
+    } // rof
+  }
+  else
+  {
+    typedef itk::ConstNeighborhoodIterator< I > TNeighIt;
+    typename I::SizeType nSize;
+    nSize.Fill( 1 );
+
+    TNeighIt nIt( nSize, this->GetInput( ), reg );
+    nIt.SetLocation( n.Vertex );
+    for( unsigned int i = 0; i < nIt.Size( ); i++ )
+    {
+      TVertex idxN = nIt.GetIndex( i );
+      if( idxN == n.Vertex )
+        continue;
+      if( !reg.IsInside( idxN ) )
+        continue;
+      N.push_back( _TNode( idxN, n.FrontId ) );
+
+    } // rof
+
+  } // fi
+}
+
+// -------------------------------------------------------------------------
+template< class I, class A >
+void fpa::Image::Algorithm< I, A >::
+_NeighsInDim( const _TNode& n, const unsigned int& d, _TNodes& N )
+{
+  typename I::RegionType reg = this->GetInput( )->GetRequestedRegion( );
+
+  N.clear( );
+  for( int i = -1; i <= 1; i += 2 )
+  {
+    TVertex v = n.Vertex;
+    v[ d ] += i;
+    if( reg.IsInside( v ) )
+      N.push_back( _TNode( v, n.FrontId ) );
+
+  } // rof
+}
+
+// -------------------------------------------------------------------------
+template< class I, class A >
+void fpa::Image::Algorithm< I, A >::
+_InitializeResults( )
+{
+}
+
+#endif // __FPA__IMAGE__ALGORITHM__HXX__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/Dijkstra.h b/lib/fpa/Image/Dijkstra.h
new file mode 100644 (file)
index 0000000..1464d3c
--- /dev/null
@@ -0,0 +1,58 @@
+#ifndef __FPA__IMAGE__DIJKSTRA__H__
+#define __FPA__IMAGE__DIJKSTRA__H__
+
+#include <itkImage.h>
+#include <itkImageFunction.h>
+#include <itkImageToImageFilter.h>
+#include <itkIndex.h>
+#include <fpa/Base/Dijkstra.h>
+#include <fpa/Image/Algorithm.h>
+
+namespace fpa
+{
+  namespace Image
+  {
+    /**
+     * @param I Input image type
+     */
+    template< class I, class C >
+    class Dijkstra
+      : public Algorithm< I, fpa::Base::Dijkstra< typename I::IndexType, C, typename I::PixelType, itk::Functor::IndexLexicographicCompare< I::ImageDimension >, itk::ImageToImageFilter< I, itk::Image< C, I::ImageDimension > > > >
+    {
+    public:
+      // Standard class typdedefs
+      typedef typename I::IndexType TVertex;
+      typedef typename I::PixelType TVertexValue;
+      typedef itk::Image< C, I::ImageDimension > TCostImage;
+      typedef itk::ImageToImageFilter< I, TCostImage > TBaseFilter;
+      typedef fpa::Base::Dijkstra< TVertex, C, TVertexValue, itk::Functor::IndexLexicographicCompare< I::ImageDimension >, TBaseFilter > TBaseAlgorithm;
+
+      typedef Dijkstra                        Self;
+      typedef Algorithm< I, TBaseAlgorithm >  Superclass;
+      typedef itk::SmartPointer< Self >       Pointer;
+      typedef itk::SmartPointer< const Self > ConstPointer;
+
+    public:
+      itkNewMacro( Self );
+      itkTypeMacro( Dijkstra, fpaBaseDijkstra );
+
+    protected:
+      Dijkstra( )
+        : Superclass( )
+        { }
+      virtual ~Dijkstra( )
+        { }
+
+    private:
+      // Purposely not implemented
+      Dijkstra( const Self& );
+      void operator=( const Self& );
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#endif // __FPA__IMAGE__DIJKSTRA__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/FastMarching.h b/lib/fpa/Image/FastMarching.h
new file mode 100644 (file)
index 0000000..d86d99d
--- /dev/null
@@ -0,0 +1,58 @@
+#ifndef __FPA__IMAGE__FASTMARCHING__H__
+#define __FPA__IMAGE__FASTMARCHING__H__
+
+#include <itkImage.h>
+#include <itkImageFunction.h>
+#include <itkImageToImageFilter.h>
+#include <itkIndex.h>
+#include <fpa/Base/FastMarching.h>
+#include <fpa/Image/Algorithm.h>
+
+namespace fpa
+{
+  namespace Image
+  {
+    /**
+     * @param I Input image type
+     */
+    template< class I, class C >
+    class FastMarching
+      : public Algorithm< I, fpa::Base::FastMarching< typename I::IndexType, C, typename I::PixelType, itk::Functor::IndexLexicographicCompare< I::ImageDimension >, itk::ImageToImageFilter< I, itk::Image< C, I::ImageDimension > > > >
+    {
+    public:
+      // Standard class typdedefs
+      typedef typename I::IndexType TVertex;
+      typedef typename I::PixelType TVertexValue;
+      typedef itk::Image< C, I::ImageDimension > TCostImage;
+      typedef itk::ImageToImageFilter< I, TCostImage > TBaseFilter;
+      typedef fpa::Base::FastMarching< TVertex, C, TVertexValue, itk::Functor::IndexLexicographicCompare< I::ImageDimension >, TBaseFilter > TBaseAlgorithm;
+
+      typedef FastMarching                        Self;
+      typedef Algorithm< I, TBaseAlgorithm >  Superclass;
+      typedef itk::SmartPointer< Self >       Pointer;
+      typedef itk::SmartPointer< const Self > ConstPointer;
+
+    public:
+      itkNewMacro( Self );
+      itkTypeMacro( FastMarching, fpaBaseFastMarching );
+
+    protected:
+      FastMarching( )
+        : Superclass( )
+        { }
+      virtual ~FastMarching( )
+        { }
+
+    private:
+      // Purposely not implemented
+      FastMarching( const Self& );
+      void operator=( const Self& );
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#endif // __FPA__IMAGE__FASTMARCHING__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/Functors/RegionGrowAllBelongsFunction.h b/lib/fpa/Image/Functors/RegionGrowAllBelongsFunction.h
new file mode 100644 (file)
index 0000000..ca52050
--- /dev/null
@@ -0,0 +1,78 @@
+#ifndef __FPA__IMAGE__FUNCTORS__REGIONGROWALLBELONGSFUNCTION__H__
+#define __FPA__IMAGE__FUNCTORS__REGIONGROWALLBELONGSFUNCTION__H__
+
+#include <itkImageFunction.h>
+
+namespace fpa
+{
+  namespace Image
+  {
+    namespace Functors
+    {
+      /**
+       */
+      template< class I >
+      class RegionGrowAllBelongsFunction
+        : public itk::ImageFunction< I, bool >
+      {
+      public:
+        /// Type-related and pointers
+        typedef RegionGrowAllBelongsFunction    Self;
+        typedef itk::ImageFunction< I, bool >   Superclass;
+        typedef itk::SmartPointer< Self >       Pointer;
+        typedef itk::SmartPointer< const Self > ConstPointer;
+
+        // Superclass' types
+        typedef typename Superclass::InputImageType      InputImageType;
+        typedef typename Superclass::InputPixelType      InputPixelType;
+        typedef typename Superclass::OutputType          OutputType;
+        typedef typename Superclass::CoordRepType        CoordRepType;
+        typedef typename Superclass::IndexType           IndexType;
+        typedef typename Superclass::IndexValueType      IndexValueType;
+        typedef typename Superclass::ContinuousIndexType ContinuousIndexType;
+        typedef typename Superclass::PointType           PointType;
+
+      public:
+        itkNewMacro( Self );
+        itkTypeMacro( RegionGrowAllBelongsFunction, itkImageFunction );
+
+      public:
+        virtual OutputType Evaluate( const PointType& point ) const
+          {
+            IndexType index;
+            this->ConvertPointToNearestIndex( point, index );
+            return( this->EvaluateAtIndex( index ) );
+          }
+        virtual OutputType EvaluateAtContinuousIndex(
+          const ContinuousIndexType& cindex
+          ) const
+          {
+            IndexType index;
+            this->ConvertContinuousIndexToNearestIndex( cindex, index );
+            return( this->EvaluateAtIndex( index ) );
+          }
+        virtual OutputType EvaluateAtIndex( const IndexType& index ) const
+          { return( true ); }
+
+      protected:
+        RegionGrowAllBelongsFunction( )
+          : Superclass( )
+          { }
+        virtual ~RegionGrowAllBelongsFunction( )
+          { }
+
+      private:
+        // Purposely not implemented
+        RegionGrowAllBelongsFunction( const Self& );
+        void operator=( const Self& );
+      };
+
+    } // ecapseman
+
+  } // ecapseman
+
+} // ecapseman
+
+#endif // __FPA__IMAGE__FUNCTORS__REGIONGROWALLBELONGSFUNCTION__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/Functors/RegionGrowThresholdFunction.h b/lib/fpa/Image/Functors/RegionGrowThresholdFunction.h
new file mode 100644 (file)
index 0000000..00a813f
--- /dev/null
@@ -0,0 +1,96 @@
+#ifndef __FPA__IMAGE__FUNCTORS__REGIONGROWTHRESHOLDFUNCTION__H__
+#define __FPA__IMAGE__FUNCTORS__REGIONGROWTHRESHOLDFUNCTION__H__
+
+#include <limits>
+#include <fpa/Image/Functors/RegionGrowAllBelongsFunction.h>
+
+namespace fpa
+{
+  namespace Image
+  {
+    namespace Functors
+    {
+      /**
+       */
+      template< class I >
+      class RegionGrowThresholdFunction
+        : public RegionGrowAllBelongsFunction< I >
+      {
+      public:
+        /// Type-related and pointers
+        typedef RegionGrowThresholdFunction       Self;
+        typedef RegionGrowAllBelongsFunction< I > Superclass;
+        typedef itk::SmartPointer< Self >         Pointer;
+        typedef itk::SmartPointer< const Self >   ConstPointer;
+
+        // Superclass' types
+        typedef typename Superclass::InputImageType      InputImageType;
+        typedef typename Superclass::InputPixelType      InputPixelType;
+        typedef typename Superclass::OutputType          OutputType;
+        typedef typename Superclass::CoordRepType        CoordRepType;
+        typedef typename Superclass::IndexType           IndexType;
+        typedef typename Superclass::IndexValueType      IndexValueType;
+        typedef typename Superclass::ContinuousIndexType ContinuousIndexType;
+        typedef typename Superclass::PointType           PointType;
+
+      public:
+        itkNewMacro( Self );
+        itkTypeMacro(
+          RegionGrowThresholdFunction,
+          RegionGrowAllBelongsFunction
+          );
+
+        itkGetConstMacro( LowerThreshold, InputPixelType );
+        itkGetConstMacro( UpperThreshold, InputPixelType );
+
+        itkSetMacro( LowerThreshold, InputPixelType );
+        itkSetMacro( UpperThreshold, InputPixelType );
+
+      public:
+        virtual OutputType EvaluateAtIndex( const IndexType& index ) const
+          {
+            const I* img = this->GetInputImage( );
+            if( img != NULL )
+            {
+              if( this->IsInsideBuffer( index ) )
+              {
+                InputPixelType v = img->GetPixel( index );
+                return(
+                  this->m_LowerThreshold <= v &&
+                  v <= this->m_UpperThreshold
+                  );
+
+              } // fi
+
+            } // fi
+            return( false );
+          }
+
+      protected:
+        RegionGrowThresholdFunction( )
+          : Superclass( ),
+            m_LowerThreshold( std::numeric_limits< InputPixelType >::min( ) ),
+            m_UpperThreshold( std::numeric_limits< InputPixelType >::max( ) )
+          { }
+        virtual ~RegionGrowThresholdFunction( )
+          { }
+
+      private:
+        // Purposely not implemented
+        RegionGrowThresholdFunction( const Self& );
+        void operator=( const Self& );
+
+      protected:
+        InputPixelType m_LowerThreshold;
+        InputPixelType m_UpperThreshold;
+      };
+
+    } // ecapseman
+
+  } // ecapseman
+
+} // ecapseman
+
+#endif // __FPA__IMAGE__FUNCTORS__REGIONGROWTHRESHOLDFUNCTION__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/RegionGrow.h b/lib/fpa/Image/RegionGrow.h
new file mode 100644 (file)
index 0000000..d65fa6c
--- /dev/null
@@ -0,0 +1,83 @@
+#ifndef __FPA__IMAGE__REGIONGROW__H__
+#define __FPA__IMAGE__REGIONGROW__H__
+
+#include <itkImageFunction.h>
+#include <itkImageToImageFilter.h>
+#include <itkIndex.h>
+#include <fpa/Base/RegionGrow.h>
+#include <fpa/Image/Algorithm.h>
+
+namespace fpa
+{
+  namespace Image
+  {
+    /**
+     * @param I Input image type
+     */
+    template< class I >
+    class RegionGrow
+      : public Algorithm< I, fpa::Base::RegionGrow< typename I::IndexType, typename I::PixelType, typename I::PixelType, itk::Functor::IndexLexicographicCompare< I::ImageDimension >, itk::ImageToImageFilter< I, I > > >
+    {
+    public:
+      // Standard class typdedefs
+      typedef typename I::IndexType TVertex;
+      typedef typename I::PixelType TResult;
+      typedef typename I::PixelType TVertexValue;
+      typedef itk::ImageToImageFilter< I, I > TBaseFilter;
+      typedef fpa::Base::RegionGrow< TVertex, TResult, TVertexValue, itk::Functor::IndexLexicographicCompare< I::ImageDimension >, TBaseFilter > TBaseAlgorithm;
+
+      typedef RegionGrow                      Self;
+      typedef Algorithm< I, TBaseAlgorithm >  Superclass;
+      typedef itk::SmartPointer< Self >       Pointer;
+      typedef itk::SmartPointer< const Self > ConstPointer;
+
+      typedef itk::ImageFunction< I, bool > TMembershipFunction;
+
+    public:
+      itkNewMacro( Self );
+      itkTypeMacro( RegionGrow, fpaBaseRegionGrow );
+
+      itkGetObjectMacro( MembershipFunction, TMembershipFunction );
+      itkSetObjectMacro( MembershipFunction, TMembershipFunction );
+
+    protected:
+      RegionGrow( )
+        : Superclass( ),
+          m_MembershipFunction( NULL )
+        { }
+      virtual ~RegionGrow( )
+        { }
+
+      virtual bool _CheckMembership(
+        const typename Superclass::_TNode& n
+        ) const
+        {
+          if( this->m_MembershipFunction.IsNotNull( ) )
+          {
+            if(
+              this->m_MembershipFunction->GetInputImage( ) !=
+              this->GetInput( )
+              )
+              this->m_MembershipFunction->SetInputImage( this->GetInput( ) );
+            return( this->m_MembershipFunction->EvaluateAtIndex( n.Vertex ) );
+          }
+          else
+            return( false );
+        }
+
+    private:
+      // Purposely not implemented
+      RegionGrow( const Self& );
+      void operator=( const Self& );
+
+    protected:
+      typename TMembershipFunction::Pointer m_MembershipFunction;
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#endif // __FPA__IMAGE__REGIONGROW__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/RegionGrowWithMultipleThresholds.h b/lib/fpa/Image/RegionGrowWithMultipleThresholds.h
new file mode 100644 (file)
index 0000000..11066fe
--- /dev/null
@@ -0,0 +1,70 @@
+#ifndef __FPA__IMAGE__REGIONGROWWITHMULTIPLETHRESHOLDS__H__
+#define __FPA__IMAGE__REGIONGROWWITHMULTIPLETHRESHOLDS__H__
+
+#include <map>
+#include <fpa/Image/RegionGrow.h>
+
+namespace fpa
+{
+  namespace Image
+  {
+    /**
+     * @param I Input image type
+     */
+    template< class I >
+    class RegionGrowWithMultipleThresholds
+      : public RegionGrow< I >
+    {
+    public:
+      typedef RegionGrowWithMultipleThresholds Self;
+      typedef RegionGrow< I >                  Superclass;
+      typedef itk::SmartPointer< Self >        Pointer;
+      typedef itk::SmartPointer< const Self >  ConstPointer;
+
+      typedef typename I::PixelType TPixel;
+
+      typedef std::map< TPixel, unsigned long > THistogram;
+      typedef typename Superclass::TBaseAlgorithm TBaseAlgorithm;
+
+    protected:
+      typedef typename TBaseAlgorithm::_TNode _TNode;
+
+    public:
+      itkNewMacro( Self );
+      itkTypeMacro( RegionGrowWithMultipleThresholds, RegionGrow );
+
+      itkGetConstMacro( DerivativeThreshold, double );
+      itkSetMacro( DerivativeThreshold, double );
+
+    public:
+      void AddThreshold( const TPixel& v );
+      void AddThresholds(
+        const TPixel& t0, const TPixel& t1,
+        const unsigned int& samples
+        );
+
+    protected:
+      RegionGrowWithMultipleThresholds( );
+      virtual ~RegionGrowWithMultipleThresholds( );
+
+      virtual bool _UpdateResult( _TNode& n );
+      virtual void _AfterLoop( );
+
+    private:
+      RegionGrowWithMultipleThresholds( const Self& ); // Not impl.
+      void operator=( const Self& );                   // Not impl.
+
+    protected:
+      double m_DerivativeThreshold;
+      THistogram m_Histogram;
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#include <fpa/Image/RegionGrowWithMultipleThresholds.hxx>
+
+#endif // __FPA__IMAGE__REGIONGROWWITHMULTIPLETHRESHOLDS__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/Image/RegionGrowWithMultipleThresholds.hxx b/lib/fpa/Image/RegionGrowWithMultipleThresholds.hxx
new file mode 100644 (file)
index 0000000..b4d4331
--- /dev/null
@@ -0,0 +1,137 @@
+#ifndef __FPA__IMAGE__REGIONGROWWITHMULTIPLETHRESHOLDS__HXX__
+#define __FPA__IMAGE__REGIONGROWWITHMULTIPLETHRESHOLDS__HXX__
+
+#include <fpa/Image/Functors/RegionGrowThresholdFunction.h>
+
+// -------------------------------------------------------------------------
+template< class I >
+void fpa::Image::RegionGrowWithMultipleThresholds< I >::
+AddThreshold( const TPixel& v )
+{
+  this->m_Histogram[ v ] = 0;
+
+  typedef
+    fpa::Image::Functors::RegionGrowThresholdFunction< I >
+    TFunction;
+  TFunction* function =
+    dynamic_cast< TFunction* >( this->GetMembershipFunction( ) );
+  if( function != NULL )
+  {
+    function->SetLowerThreshold( this->m_Histogram.begin( )->first );
+    function->SetUpperThreshold( this->m_Histogram.rbegin( )->first );
+
+  } // fi
+  this->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class I >
+void fpa::Image::RegionGrowWithMultipleThresholds< I >::
+AddThresholds(
+  const TPixel& t0,
+  const TPixel& t1,
+  const unsigned int& samples
+  )
+{
+  double off =
+    ( double( t1 ) - double( t0 ) ) / ( double( samples ) - double( 1 ) );
+  for( unsigned int s = 0; s < samples; s++ )
+    this->AddThreshold( TPixel( ( double( s ) * off ) + double( t0 ) ) );
+}
+
+// -------------------------------------------------------------------------
+template< class I >
+fpa::Image::RegionGrowWithMultipleThresholds< I >::
+RegionGrowWithMultipleThresholds( )
+  : Superclass( ),
+    m_DerivativeThreshold( double( 3 ) )
+{
+  typedef
+    fpa::Image::Functors::RegionGrowThresholdFunction< I >
+    TFunction;
+  typename TFunction::Pointer function = TFunction::New( );
+  this->SetMembershipFunction( function );
+}
+
+// -------------------------------------------------------------------------
+template< class I >
+fpa::Image::RegionGrowWithMultipleThresholds< I >::
+~RegionGrowWithMultipleThresholds( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class I >
+bool fpa::Image::RegionGrowWithMultipleThresholds< I >::
+_UpdateResult( _TNode& n )
+{
+  bool ret = this->Superclass::_UpdateResult( n );
+
+  if( ret )
+  {
+    TPixel v = TPixel( this->_Cost( n.Vertex, n.Vertex ) );
+
+    typename THistogram::reverse_iterator hIt = this->m_Histogram.rbegin( );
+    while( hIt != this->m_Histogram.rend( ) )
+    {
+      if( v <= hIt->first )
+      {
+        hIt->second += 1;
+        hIt++;
+      }
+      else
+        hIt = this->m_Histogram.rend( );
+
+    } // elihw
+    this->GetOutput( )->SetPixel( n.Vertex, v );
+
+  } // fi
+  return( ret );
+}
+
+// -------------------------------------------------------------------------
+template< class I >
+void fpa::Image::RegionGrowWithMultipleThresholds< I >::
+_AfterLoop( )
+{
+  typename THistogram::iterator prevIt = this->m_Histogram.begin( );
+  typename THistogram::iterator currIt = prevIt; currIt++;
+  typename THistogram::iterator nextIt = currIt; nextIt++;
+  double prev_d1 = double( 0 );
+  for( ; nextIt != this->m_Histogram.end( ); ++prevIt, ++currIt, ++nextIt )
+  {
+    double d1 = double( nextIt->second ) - double( prevIt->second );
+    d1 /= double( nextIt->first ) - double( prevIt->first );
+
+    std::cout
+      << currIt->first  << " "
+      << currIt->second << " "
+      << d1 << " "
+      << ( d1 - prev_d1 ) << std::endl;
+
+    prev_d1 = d1;
+
+  } // rof
+
+
+  /*
+    double prev = double( 0 );
+    while( hIt != this->m_Histogram.end( ) )
+    {
+    double curr = double( hIt->second );
+    double deri = curr - prev;
+
+    if( hIt != this->m_Histogram.begin( ) )
+    std::cout << hIt->first << " " << curr << " " << deri << std::endl;
+
+    prev = curr;
+    hIt++;
+
+    } // rof
+  */
+  this->Superclass::_AfterLoop( );
+}
+
+#endif // __FPA__IMAGE__REGIONGROWWITHMULTIPLETHRESHOLDS__HXX__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/VTK/Image2DObserver.h b/lib/fpa/VTK/Image2DObserver.h
new file mode 100644 (file)
index 0000000..68dcae0
--- /dev/null
@@ -0,0 +1,91 @@
+#ifndef __FPA__VTK__IMAGE2DOBSERVER__H__
+#define __FPA__VTK__IMAGE2DOBSERVER__H__
+
+#include <set>
+#include <itkCommand.h>
+#include <vtkActor.h>
+#include <vtkImageActor.h>
+#include <vtkImageData.h>
+#include <vtkPolyData.h>
+#include <vtkSmartPointer.h>
+
+namespace fpa
+{
+  namespace VTK
+  {
+    /**
+     */
+    template< class F, class R >
+    class Image2DObserver
+      : public itk::Command
+    {
+    public:
+      typedef Image2DObserver                   Self;
+      typedef itk::Command                    Superclass;
+      typedef itk::SmartPointer< Self >       Pointer;
+      typedef itk::SmartPointer< const Self > ConstPointer;
+
+      typedef F                             TFilter;
+      typedef R                             TRenderWindow;
+      typedef typename TFilter::TInputImage TImage;
+      typedef typename TFilter::TVertex     TVertex;
+
+      typedef std::set< TVertex > TVertices;
+
+    public:
+      itkNewMacro( Self );
+      itkTypeMacro( Image2DObserver, itkCommand );
+
+    public:
+      void SetImage( const TImage* img, R* rw );
+      void SetPixel(
+        typename TImage::IndexType idx,
+        unsigned char red,
+        unsigned char green,
+        unsigned char blue,
+        unsigned char alpha
+        );
+      void Render( );
+      void Execute( itk::Object* c, const itk::EventObject& e )
+        { this->Execute( ( const itk::Object* )( c ), e ); }
+      void Execute( const itk::Object* c, const itk::EventObject& e );
+
+    protected:
+      Image2DObserver( )
+        : Superclass( ),
+          m_RenderWindow( NULL ),
+          m_Number( 0 ),
+          m_Percentage( 0 )
+        {
+          this->m_Stencil = TImageData::New( );
+          this->m_StencilActor = TImageActor::New( );
+        }
+      virtual ~Image2DObserver( )
+        { }
+
+    private:
+      Image2DObserver( const Self& ); // Not impl.
+      void operator=( const Self& );  // Not impl.
+
+    protected:
+      typedef vtkSmartPointer< vtkImageActor > TImageActor;
+      typedef vtkSmartPointer< vtkImageData >  TImageData;
+
+      typename TImage::ConstPointer m_Image;
+
+      R*            m_RenderWindow;
+      TImageData    m_Stencil;
+      TImageActor   m_StencilActor;
+      unsigned long m_Number;
+      unsigned long m_Percentage;
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#include <fpa/VTK/Image2DObserver.hxx>
+
+#endif // __FPA__VTK__IMAGE2DOBSERVER__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/VTK/Image2DObserver.hxx b/lib/fpa/VTK/Image2DObserver.hxx
new file mode 100644 (file)
index 0000000..7ebe526
--- /dev/null
@@ -0,0 +1,188 @@
+#ifndef __FPA__VTK__IMAGE2DOBSERVER__HXX__
+#define __FPA__VTK__IMAGE2DOBSERVER__HXX__
+
+#include <vtkCellArray.h>
+#include <vtkPolyDataMapper.h>
+#include <vtkPoints.h>
+#include <vtkPointData.h>
+#include <vtkProperty.h>
+#include <vtkRenderer.h>
+#include <vtkRendererCollection.h>
+
+// -------------------------------------------------------------------------
+template< class F, class R >
+void fpa::VTK::Image2DObserver< F, R >::
+SetImage( const TImage* img, R* rw )
+{
+  this->m_Image = img;
+  this->m_RenderWindow = rw;
+
+  unsigned int minD = TImage::ImageDimension;
+  minD = ( minD < 3 )? minD: 3;
+
+  int e[ 6 ] = { 0 };
+  typename TImage::RegionType reg = this->m_Image->GetRequestedRegion( );
+  for( unsigned int i = 0; i < minD; i++ )
+  {
+    e[ ( i << 1 ) + 0 ] = reg.GetIndex( )[ i ];
+    e[ ( i << 1 ) + 1 ] = reg.GetIndex( )[ i ] + reg.GetSize( )[ i ] - 1;
+
+  } // rof
+
+  typename TImage::SpacingType spac = this->m_Image->GetSpacing( );
+  double s[ 3 ] = { 1, 1, 1 };
+  for( unsigned int i = 0; i < minD; i++ )
+    s[ i ] = double( spac[ i ] );
+
+  typename TImage::PointType orig = this->m_Image->GetOrigin( );
+  double o[ 3 ] = { 0 };
+  for( unsigned int i = 0; i < minD; i++ )
+    o[ i ] = double( orig[ i ] );
+
+  this->m_Stencil->SetExtent( e );
+  this->m_Stencil->SetSpacing( s );
+  this->m_Stencil->SetOrigin( o );
+  this->m_Stencil->AllocateScalars( VTK_UNSIGNED_CHAR, 4 );
+  for( unsigned int i = 0; i < 3; i++ )
+    this->m_Stencil->GetPointData( )->
+      GetScalars( )->FillComponent( i, 255 );
+  this->m_Stencil->GetPointData( )->GetScalars( )->FillComponent( 3, 0 );
+
+  this->m_StencilActor->SetInputData( this->m_Stencil );
+  this->m_StencilActor->InterpolateOff( );
+
+  double nPix =
+    double( this->m_Image->GetRequestedRegion( ).GetNumberOfPixels( ) );
+  this->m_Percentage = ( unsigned long )( nPix * 0.01 );
+  this->m_Number = 0;
+
+  if( this->m_RenderWindow != NULL )
+  {
+    vtkRenderer* ren =
+      this->m_RenderWindow->GetRenderers( )->GetFirstRenderer( );
+    ren->AddActor( this->m_StencilActor );
+    this->Render( );
+
+  } // fi
+}
+
+// -------------------------------------------------------------------------
+template< class F, class R >
+void fpa::VTK::Image2DObserver< F, R >::
+SetPixel(
+  typename TImage::IndexType idx,
+  unsigned char red,
+  unsigned char green,
+  unsigned char blue,
+  unsigned char alpha
+  )
+{
+  this->m_Stencil->SetScalarComponentFromDouble
+    ( idx[ 0 ], idx[ 1 ], 0, 0, red );
+  this->m_Stencil->SetScalarComponentFromDouble
+    ( idx[ 0 ], idx[ 1 ], 0, 1, green );
+  this->m_Stencil->SetScalarComponentFromDouble
+    ( idx[ 0 ], idx[ 1 ], 0, 2, blue );
+  this->m_Stencil->SetScalarComponentFromDouble
+    ( idx[ 0 ], idx[ 1 ], 0, 3, alpha );
+  this->m_Stencil->Modified( );
+  this->m_StencilActor->Modified( );
+}
+
+// -------------------------------------------------------------------------
+template< class F, class R >
+void fpa::VTK::Image2DObserver< F, R >::
+Render( )
+{
+  if( this->m_RenderWindow != NULL )
+    this->m_RenderWindow->Render( );
+}
+
+// -------------------------------------------------------------------------
+template< class F, class R >
+void fpa::VTK::Image2DObserver< F, R >::
+Execute( const itk::Object* c, const itk::EventObject& e )
+{
+  typedef typename TImage::IndexType        TIndex;
+  typedef typename TImage::PointType        TPoint;
+  typedef typename TFilter::TEvent          TEvent;
+  typedef typename TFilter::TFrontEvent     TFrontEvent;
+  typedef typename TFilter::TMarkEvent      TMarkEvent;
+  typedef typename TFilter::TCollisionEvent TCollisionEvent;
+  typedef typename TFilter::TEndEvent       TEndEvent;
+
+  static unsigned char Colors[][4] =
+    {
+      {   0,   0, 127, 127 },
+      {   0, 127, 127, 127 },
+      { 127,   0, 127, 127 },
+      { 127, 127,   0, 127 },
+      {   0,   0,  63, 127 },
+      {   0,  63,  63, 127 },
+      {  63,   0,  63, 127 },
+      {  63,  63,   0, 127 },
+      {  63,  63, 127, 127 },
+      {  63, 127, 127, 127 },
+      { 127,  63, 127, 127 },
+      { 127, 127,  63, 127 },
+      { 127, 127,  63, 127 },
+      { 127,  63,  63, 127 },
+      {  63, 127,  63, 127 },
+      {  63,  63, 127, 127 }
+    };
+
+  if( this->m_RenderWindow == NULL )
+    return;
+
+  const TEvent* evt = dynamic_cast< const TEvent* >( &e );
+  TFrontEvent fevt;
+  TMarkEvent mevt;
+  TCollisionEvent cevt;
+  TEndEvent eevt;
+  if( evt != NULL )
+  {
+    if( mevt.CheckEvent( evt ) )
+      this->SetPixel(
+        evt->Node.Vertex,
+        Colors[ evt->Node.FrontId ][ 0 ],
+        Colors[ evt->Node.FrontId ][ 1 ],
+        Colors[ evt->Node.FrontId ][ 2 ],
+        Colors[ evt->Node.FrontId ][ 3 ]
+        );
+    else if( fevt.CheckEvent( evt ) )
+      this->SetPixel( evt->Node.Vertex, 255, 0, 0, 255 );
+    if( cevt.CheckEvent( evt ) )
+      this->SetPixel( evt->Node.Vertex, 100, 50, 25, 255 );
+
+    // Update visualization
+    if( this->m_Number == 0 || this->m_Percentage < 0 )
+      this->Render( );
+    this->m_Number++;
+    this->m_Number %= this->m_Percentage;
+
+    if( eevt.CheckEvent( evt ) )
+    {
+      if( this->m_RenderWindow != NULL )
+      {
+        /* TODO
+           vtkRenderer* ren =
+           this->m_RenderWindow->GetRenderers( )->GetFirstRenderer( );
+           ren->RemoveActor( this->m_StencilActor );
+        */
+        this->Render( );
+
+      } // fi
+
+    } // fi
+  }
+  else
+  {
+    // TODO: Remove all!
+    this->Render( );
+
+  } // fi
+}
+
+#endif // __FPA__VTK__IMAGE2DOBSERVER__HXX__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/VTK/Image3DObserver.h b/lib/fpa/VTK/Image3DObserver.h
new file mode 100644 (file)
index 0000000..860946d
--- /dev/null
@@ -0,0 +1,93 @@
+#ifndef __FPA__VTK__IMAGE3DOBSERVER__H__
+#define __FPA__VTK__IMAGE3DOBSERVER__H__
+
+#include <map>
+
+#include <itkCommand.h>
+
+#include <vtkActor.h>
+#include <vtkPolyData.h>
+#include <vtkPolyDataMapper.h>
+#include <vtkSmartPointer.h>
+
+namespace fpa
+{
+  namespace VTK
+  {
+    /**
+     */
+    template< class F, class R >
+    class Image3DObserver
+      : public itk::Command
+    {
+    public:
+      typedef Image3DObserver                 Self;
+      typedef itk::Command                    Superclass;
+      typedef itk::SmartPointer< Self >       Pointer;
+      typedef itk::SmartPointer< const Self > ConstPointer;
+
+      typedef F                             TFilter;
+      typedef R                             TRenderWindow;
+      typedef typename TFilter::TInputImage TImage;
+      typedef typename TFilter::TVertex     TVertex;
+
+      typedef std::pair< unsigned long, unsigned long > TVertexIds;
+
+      struct TVertexCmp
+      {
+        bool operator()( const TVertex& a, const TVertex& b ) const
+          {
+            unsigned int d = 0;
+            while( d < TVertex::Dimension && a[ d ] == b[ d ] )
+              d++;
+            if( d < TVertex::Dimension )
+              return( a[ d ] < b[ d ] );
+            else
+              return( false );
+          }
+      };
+      typedef std::map< TVertex, TVertexIds, TVertexCmp > TVertices;
+
+    public:
+      itkNewMacro( Self );
+      itkTypeMacro( Image3DObserver, itkCommand );
+
+    public:
+      void SetImage( const TImage* img, R* rw );
+      void Render( );
+      void Execute( itk::Object* c, const itk::EventObject& e )
+        { this->Execute( ( const itk::Object* )( c ), e ); }
+      void Execute( const itk::Object* c, const itk::EventObject& e );
+
+    protected:
+      Image3DObserver( );
+      virtual ~Image3DObserver( )
+        { }
+
+    private:
+      Image3DObserver( const Self& ); // Not impl.
+      void operator=( const Self& );  // Not impl.
+
+    protected:
+      typename TImage::ConstPointer m_Image;
+
+      R*            m_RenderWindow;
+      unsigned long m_Number;
+      unsigned long m_Percentage;
+
+      vtkSmartPointer< vtkPolyData > m_PolyData;
+      vtkSmartPointer< vtkPolyDataMapper > m_PolyDataMapper;
+      vtkSmartPointer< vtkActor > m_PolyDataActor;
+
+      TVertices m_Vertices;
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#include <fpa/VTK/Image3DObserver.hxx>
+
+#endif // __FPA__VTK__IMAGE3DOBSERVER__H__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/VTK/Image3DObserver.hxx b/lib/fpa/VTK/Image3DObserver.hxx
new file mode 100644 (file)
index 0000000..e79bbbe
--- /dev/null
@@ -0,0 +1,203 @@
+#ifndef __FPA__VTK__IMAGE3DOBSERVER__HXX__
+#define __FPA__VTK__IMAGE3DOBSERVER__HXX__
+
+#include <vtkCellArray.h>
+#include <vtkPoints.h>
+#include <vtkProperty.h>
+#include <vtkRenderer.h>
+#include <vtkRendererCollection.h>
+/*
+  #include <vtkCellArray.h>
+  #include <vtkPolyDataMapper.h>
+  #include <vtkPoints.h>
+  #include <vtkPointData.h>
+*/
+
+// -------------------------------------------------------------------------
+template< class F, class R >
+void fpa::VTK::Image3DObserver< F, R >::
+SetImage( const TImage* img, R* rw )
+{
+  this->m_Image = img;
+  this->m_RenderWindow = rw;
+  this->m_Vertices.clear( );
+
+  /*
+    unsigned int minD = TImage::ImageDimension;
+    minD = ( minD < 3 )? minD: 3;
+
+    int e[ 6 ] = { 0 };
+    typename TImage::RegionType reg = this->m_Image->GetRequestedRegion( );
+    for( unsigned int i = 0; i < minD; i++ )
+    {
+    e[ ( i << 1 ) + 0 ] = reg.GetIndex( )[ i ];
+    e[ ( i << 1 ) + 1 ] = reg.GetIndex( )[ i ] + reg.GetSize( )[ i ] - 1;
+
+    } // rof
+
+    typename TImage::SpacingType spac = this->m_Image->GetSpacing( );
+    double s[ 3 ] = { 1, 1, 1 };
+    for( unsigned int i = 0; i < minD; i++ )
+    s[ i ] = double( spac[ i ] );
+
+    typename TImage::PointType orig = this->m_Image->GetOrigin( );
+    double o[ 3 ] = { 0 };
+    for( unsigned int i = 0; i < minD; i++ )
+    o[ i ] = double( orig[ i ] );
+
+    this->m_Stencil->SetExtent( e );
+    this->m_Stencil->SetSpacing( s );
+    this->m_Stencil->SetOrigin( o );
+    this->m_Stencil->AllocateScalars( VTK_UNSIGNED_CHAR, 4 );
+    for( unsigned int i = 0; i < 3; i++ )
+    this->m_Stencil->GetPointData( )->
+    GetScalars( )->FillComponent( i, 255 );
+    this->m_Stencil->GetPointData( )->GetScalars( )->FillComponent( 3, 0 );
+
+    this->m_StencilActor->SetInputData( this->m_Stencil );
+    this->m_StencilActor->InterpolateOff( );
+  */
+
+  double nPix =
+    double( this->m_Image->GetRequestedRegion( ).GetNumberOfPixels( ) );
+  this->m_Percentage = ( unsigned long )( nPix * 0.001 );
+  this->m_Number = 0;
+
+  if( this->m_RenderWindow != NULL )
+  {
+    vtkRenderer* ren =
+      this->m_RenderWindow->GetRenderers( )->GetFirstRenderer( );
+    ren->AddActor( this->m_PolyDataActor );
+    this->Render( );
+
+  } // fi
+}
+
+// -------------------------------------------------------------------------
+template< class F, class R >
+void fpa::VTK::Image3DObserver< F, R >::
+Render( )
+{
+  if( this->m_RenderWindow != NULL )
+    this->m_RenderWindow->Render( );
+}
+
+// -------------------------------------------------------------------------
+template< class F, class R >
+void fpa::VTK::Image3DObserver< F, R >::
+Execute( const itk::Object* c, const itk::EventObject& e )
+{
+  typedef typename TImage::IndexType        TIndex;
+  typedef typename TImage::PointType        TPoint;
+  typedef typename TFilter::TEvent          TEvent;
+  typedef typename TFilter::TFrontEvent     TFrontEvent;
+  typedef typename TFilter::TMarkEvent      TMarkEvent;
+  typedef typename TFilter::TCollisionEvent TCollisionEvent;
+  typedef typename TFilter::TEndEvent       TEndEvent;
+
+  if( this->m_RenderWindow == NULL )
+    return;
+
+  const TEvent* evt = dynamic_cast< const TEvent* >( &e );
+  TFrontEvent fevt;
+  TMarkEvent mevt;
+  TCollisionEvent cevt;
+  TEndEvent eevt;
+  if( evt != NULL )
+  {
+    typename TImage::PointType pnt;
+    this->m_Image->TransformIndexToPhysicalPoint( evt->Node.Vertex, pnt );
+
+    if( fevt.CheckEvent( evt ) )
+    {
+      if(
+        this->m_Vertices.find( evt->Node.Vertex ) == this->m_Vertices.end( )
+        )
+      {
+        unsigned long pId = this->m_PolyData->GetPoints( )->InsertNextPoint(
+          double( pnt[ 0 ] ), double( pnt[ 1 ] ), double( pnt[ 2 ] )
+          );
+        unsigned long cId = this->m_PolyData->GetVerts( )->InsertNextCell( 1 );
+        this->m_PolyData->GetVerts( )->InsertCellPoint( pId );
+        this->m_PolyData->Modified( );
+        this->m_Vertices[ evt->Node.Vertex ] = TVertexIds( pId, cId );
+
+      } // rof
+    }
+    else if( mevt.CheckEvent( evt ) )
+    {
+      /*
+        typename TVertices::iterator vIt =
+        this->m_Vertices.find( evt->Node.Vertex );
+        if( vIt != this->m_Vertices.end( ) )
+        {
+        std::cout << "Erase " << evt->Node.Vertex << " " << vIt->second.second << std::endl;
+        } // fi
+      */
+    } // fi
+
+    if( cevt.CheckEvent( evt ) )
+    {
+      // std::cout << "Collision";
+    } // fi
+    /*
+      this->SetPixel(
+      evt->Node.Vertex,
+      Colors[ evt->Node.FrontId ][ 0 ],
+      Colors[ evt->Node.FrontId ][ 1 ],
+      Colors[ evt->Node.FrontId ][ 2 ],
+      Colors[ evt->Node.FrontId ][ 3 ]
+      );
+      else if( mevt.CheckEvent( evt ) )
+      this->SetPixel( evt->Node.Vertex, 255, 0, 0, 255 );
+
+      if( cevt.CheckEvent( evt ) )
+      this->SetPixel( evt->Node.Vertex, 100, 50, 25, 255 );
+    */
+    
+    // Update visualization
+    if( this->m_Number == 0 || this->m_Percentage < 0 )
+      this->Render( );
+    this->m_Number++;
+    this->m_Number %= this->m_Percentage;
+
+    if( eevt.CheckEvent( evt ) )
+    {
+      if( this->m_RenderWindow != NULL )
+      {
+        vtkRenderer* ren =
+          this->m_RenderWindow->GetRenderers( )->GetFirstRenderer( );
+        ren->RemoveActor( this->m_PolyDataActor );
+        this->Render( );
+
+      } // fi
+
+    } // fi
+
+  } // fi
+}
+
+// -------------------------------------------------------------------------
+template< class F, class R >
+fpa::VTK::Image3DObserver< F, R >::
+Image3DObserver( )
+  : Superclass( ),
+    m_RenderWindow( NULL ),
+    m_Number( 0 ),
+    m_Percentage( 0 )
+{
+  this->m_PolyData = vtkSmartPointer< vtkPolyData >::New( );
+  this->m_PolyDataMapper = vtkSmartPointer< vtkPolyDataMapper >::New( );
+  this->m_PolyDataActor = vtkSmartPointer< vtkActor >::New( );
+
+  this->m_PolyData->SetPoints( vtkPoints::New( ) );
+  this->m_PolyData->SetVerts( vtkCellArray::New( ) );
+
+  this->m_PolyDataMapper->SetInputData( this->m_PolyData );
+  this->m_PolyDataActor->SetMapper( this->m_PolyDataMapper );
+  this->m_PolyDataActor->GetProperty( )->SetColor( 0, 0, 1 );
+}
+
+#endif // __FPA__VTK__IMAGE3DOBSERVER__HXX__
+
+// eof - $RCSfile$
diff --git a/lib/fpa/VTK/ImageMPR.cxx b/lib/fpa/VTK/ImageMPR.cxx
new file mode 100644 (file)
index 0000000..d435bb8
--- /dev/null
@@ -0,0 +1,159 @@
+#include <fpa/VTK/ImageMPR.h>
+
+#include <vtkInteractorStyleSwitch.h>
+
+// -------------------------------------------------------------------------
+fpa::VTK::ImageMPR::
+ImageMPR( )
+{
+  // vtkSmartPointer< vtkImageData > m_Image;
+
+  this->m_Outline = vtkSmartPointer< vtkOutlineSource >::New( );
+  this->m_OutlineMapper = vtkSmartPointer< vtkPolyDataMapper >::New( );
+  this->m_OutlineActor = vtkSmartPointer< vtkActor >::New( );
+  this->m_Picker = vtkSmartPointer< vtkCellPicker >::New( );
+  this->m_WidgetX = vtkSmartPointer< vtkImagePlaneWidget >::New( );
+  this->m_WidgetY = vtkSmartPointer< vtkImagePlaneWidget >::New( );
+  this->m_WidgetZ = vtkSmartPointer< vtkImagePlaneWidget >::New( );
+  this->m_Renderer = vtkSmartPointer< vtkRenderer >::New( );
+  this->m_Window = vtkSmartPointer< vtkRenderWindow >::New( );
+  this->m_Interactor = vtkSmartPointer< vtkRenderWindowInteractor >::New( );
+}
+    
+// -------------------------------------------------------------------------
+fpa::VTK::ImageMPR::
+~ImageMPR( )
+{
+}
+
+// -------------------------------------------------------------------------
+void fpa::VTK::ImageMPR::
+SetImage( vtkImageData* image )
+{
+  this->m_Image = image;
+
+  // Outline
+  this->m_Outline->SetBounds( this->m_Image->GetBounds( ) );
+  this->m_OutlineMapper->
+    SetInputConnection( this->m_Outline->GetOutputPort( ) );
+  this->m_OutlineActor->SetMapper( this->m_OutlineMapper );
+  this->m_OutlineActor->GetProperty( )->SetColor( 1, 1, 1 );
+
+  // Local picker
+  this->m_Picker->SetTolerance( 0.005 );
+
+  // Image planes
+  this->m_WidgetX->DisplayTextOn( );
+  this->m_WidgetX->SetInputData( this->m_Image );
+  this->m_WidgetX->SetPlaneOrientationToXAxes( );
+  this->m_WidgetX->SetSliceIndex( this->m_Image->GetExtent( )[ 0 ] );
+  this->m_WidgetX->SetPicker( this->m_Picker );
+  this->m_WidgetX->SetKeyPressActivationValue( 'x' );
+  this->m_WidgetX->GetPlaneProperty( )->SetLineWidth( 3 );
+  this->m_WidgetX->GetPlaneProperty( )->SetColor( 0, 1, 1 );
+
+  this->m_WidgetY->DisplayTextOn( );
+  this->m_WidgetY->SetInputData( this->m_Image );
+  this->m_WidgetY->SetPlaneOrientationToYAxes( );
+  this->m_WidgetY->SetSliceIndex( this->m_Image->GetExtent( )[ 2 ] );
+  this->m_WidgetY->SetPicker( this->m_Picker );
+  this->m_WidgetY->SetKeyPressActivationValue( 'y' );
+  this->m_WidgetY->GetPlaneProperty( )->SetColor( 1, 0, 1 );
+  this->m_WidgetY->GetPlaneProperty( )->SetLineWidth( 3 );
+  this->m_WidgetY->SetLookupTable( this->m_WidgetX->GetLookupTable( ) );
+
+  this->m_WidgetZ->DisplayTextOn( );
+  this->m_WidgetZ->SetInputData( this->m_Image );
+  this->m_WidgetZ->SetPlaneOrientationToZAxes( );
+  this->m_WidgetZ->SetSliceIndex( this->m_Image->GetExtent( )[ 4 ] );
+  this->m_WidgetZ->SetPicker( this->m_Picker );
+  this->m_WidgetZ->SetKeyPressActivationValue( 'z' );
+  this->m_WidgetZ->GetPlaneProperty( )->SetColor( 1, 1, 0 );
+  this->m_WidgetZ->GetPlaneProperty( )->SetLineWidth( 3 );
+  this->m_WidgetZ->SetLookupTable( this->m_WidgetX->GetLookupTable( ) );
+
+  // Rendering stuff
+  this->m_Window->AddRenderer( this->m_Renderer );
+  this->m_Interactor->SetRenderWindow( this->m_Window );
+
+  // Command
+  double spac[ 3 ];
+  this->m_Image->GetSpacing( spac );
+  double min_spacing = spac[ 0 ];
+  for( unsigned int d = 1; d < 3; d++ )
+    min_spacing = ( spac[ d ] < min_spacing )? spac[ d ]: min_spacing;
+
+  vtkInteractorStyleSwitch* iswitch =
+    dynamic_cast< vtkInteractorStyleSwitch* >(
+      this->m_Interactor->GetInteractorStyle( )
+      );
+  if( iswitch != NULL )
+    iswitch->SetCurrentStyleToTrackballCamera( );
+
+  // Add actors
+  this->m_Renderer->AddActor( this->m_OutlineActor );
+
+  // Prepare widgets
+  this->m_WidgetX->SetInteractor( this->m_Interactor );
+  this->m_WidgetY->SetInteractor( this->m_Interactor );
+  this->m_WidgetZ->SetInteractor( this->m_Interactor );
+}
+
+// -------------------------------------------------------------------------
+void fpa::VTK::ImageMPR::
+SetBackground( double r, double g, double b )
+{
+  this->m_Renderer->SetBackground( r, g, b );
+}
+
+// -------------------------------------------------------------------------
+void fpa::VTK::ImageMPR::
+SetSize( unsigned int w, unsigned int h )
+{
+  this->m_Window->SetSize( w, h );
+}
+
+// -------------------------------------------------------------------------
+void fpa::VTK::ImageMPR::
+AddPolyData( vtkPolyData* pd, double r, double g, double b )
+{
+  unsigned int i = this->m_PolyDatas.size( );
+
+  this->m_PolyDatas.push_back( pd );
+  this->m_Mappers.push_back( vtkSmartPointer< vtkPolyDataMapper >::New( ) );
+  this->m_Actors.push_back( vtkSmartPointer< vtkActor >::New( ) );
+
+  this->m_Mappers[ i ]->SetInputData( pd );
+  this->m_Actors[ i ]->SetMapper( this->m_Mappers[ i ] );
+  this->m_Actors[ i ]->GetProperty( )->SetColor( r, g, b );
+  this->m_Renderer->AddActor( this->m_Actors[ i ] );
+}
+
+// -------------------------------------------------------------------------
+vtkRenderWindow* fpa::VTK::ImageMPR::
+GetWindow( ) const
+{
+  return( this->m_Window );
+}
+
+// -------------------------------------------------------------------------
+vtkRenderer* fpa::VTK::ImageMPR::
+GetRenderer( ) const
+{
+  return( this->m_Renderer );
+}
+
+// -------------------------------------------------------------------------
+void fpa::VTK::ImageMPR::
+Start( )
+{
+  this->m_WidgetX->On( );
+  this->m_WidgetY->On( );
+  this->m_WidgetZ->On( );
+
+  this->m_Renderer->ResetCamera( );
+  this->m_Interactor->Initialize( );
+  this->m_Interactor->Start( );
+}
+
+// eof - $RCSfile$
diff --git a/lib/fpa/VTK/ImageMPR.h b/lib/fpa/VTK/ImageMPR.h
new file mode 100644 (file)
index 0000000..e7a65d3
--- /dev/null
@@ -0,0 +1,67 @@
+#ifndef __FPA__VTK__IMAGEMPR__H__
+#define __FPA__VTK__IMAGEMPR__H__
+
+#include <fpa/FrontAlgorithms_Export.h>
+
+#include <vtkActor.h>
+#include <vtkCellPicker.h>
+#include <vtkImageData.h>
+#include <vtkImagePlaneWidget.h>
+#include <vtkOutlineSource.h>
+#include <vtkPolyDataMapper.h>
+#include <vtkProperty.h>
+#include <vtkRenderer.h>
+#include <vtkRendererCollection.h>
+#include <vtkRenderWindow.h>
+#include <vtkRenderWindowInteractor.h>
+#include <vtkSmartPointer.h>
+#include <vtkSphereSource.h>
+
+namespace fpa
+{
+  namespace VTK
+  {
+    /**
+     */
+    class FrontAlgorithms_EXPORT ImageMPR
+    {
+    public:
+      ImageMPR( );
+      virtual ~ImageMPR( );
+
+      void SetImage( vtkImageData* image );
+      void SetBackground( double r, double g, double b );
+      void SetSize( unsigned int w, unsigned int h );
+
+      void AddPolyData( vtkPolyData* pd, double r, double g, double b );
+
+      vtkRenderWindow* GetWindow( ) const;
+      vtkRenderer* GetRenderer( ) const;
+
+      void Start( );
+
+    protected:
+      vtkSmartPointer< vtkImageData > m_Image;
+      vtkSmartPointer< vtkOutlineSource > m_Outline;
+      vtkSmartPointer< vtkPolyDataMapper > m_OutlineMapper;
+      vtkSmartPointer< vtkActor > m_OutlineActor;
+      vtkSmartPointer< vtkCellPicker > m_Picker;
+      vtkSmartPointer< vtkImagePlaneWidget > m_WidgetX;
+      vtkSmartPointer< vtkImagePlaneWidget > m_WidgetY;
+      vtkSmartPointer< vtkImagePlaneWidget > m_WidgetZ;
+      vtkSmartPointer< vtkRenderer > m_Renderer;
+      vtkSmartPointer< vtkRenderWindow > m_Window;
+      vtkSmartPointer< vtkRenderWindowInteractor > m_Interactor;
+
+      std::vector< vtkSmartPointer< vtkPolyData > > m_PolyDatas;
+      std::vector< vtkSmartPointer< vtkPolyDataMapper > > m_Mappers;
+      std::vector< vtkSmartPointer< vtkActor > > m_Actors;
+    };
+
+  } // ecapseman
+
+} // ecapseman
+
+#endif // __FPA__VTK__IMAGEMPR__H__
+
+// eof - $RCSfile$