appli/examples/example_ImageAlgorithmDijkstra_01.cxx

   1 #include <cmath>
   2 #include <iostream>
   3 #include <limits>
   4 #include <string>
   5
   6 #include <itkImage.h>
   7 #include <itkImageFileReader.h>
   8 #include <itkImageToVTKImageFilter.h>
   9 #include <itkMinimumMaximumImageCalculator.h>
  10 #include <itkInvertIntensityImageFilter.h>
  11 #include <itkDanielssonDistanceMapImageFilter.h>
  12 #include <itkUnaryFunctorImageFilter.h>
  13
  14 #include <vtkImageActor.h>
  15 #include <vtkInteractorStyleImage.h>
  16 #include <vtkPointHandleRepresentation3D.h>
  17 #include <vtkProperty.h>
  18 #include <vtkRenderer.h>
  19 #include <vtkRenderWindow.h>
  20 #include <vtkRenderWindowInteractor.h>
  21 #include <vtkSeedRepresentation.h>
  22 #include <vtkSeedWidget.h>
  23 #include <vtkSmartPointer.h>
  24
  25 #include <fpa/Image/Dijkstra.h>
  26 #include <fpa/VTK/Image2DObserver.h>
  27
  28 // -------------------------------------------------------------------------
  29 const unsigned int Dim = 2;
  30 typedef unsigned char TPixel;
  31 typedef double TScalar;
  32 typedef itk::Image< TPixel, Dim >  TImage;
  33 typedef itk::Image< TScalar, Dim >  TScalarImage;
  34 typedef itk::ImageToVTKImageFilter< TImage > TVTKImage;
  35
  36 typedef itk::ImageFileReader< TImage > TImageReader;
  37 typedef fpa::Image::Dijkstra< TScalarImage, TScalar > TFrontAlgorithm;
  38
  39 typedef
  40 fpa::VTK::Image2DObserver< TFrontAlgorithm, vtkRenderWindow >
  41 TObserver;
  42
  43 // -------------------------------------------------------------------------
  44 class InvertPixelFunctor
  45 {
  46 public:
  47   InvertPixelFunctor( )
  48     { }
  49   virtual ~InvertPixelFunctor( )
  50     { }
  51   bool operator!=( const InvertPixelFunctor& other ) const
  52     { return( false ); }
  53   bool operator==( const InvertPixelFunctor& other ) const
  54     { return( !( *this != other ) ); }
  55   inline TScalar operator()( const TScalar& A ) const
  56     {
  57       return( TScalar( 1 ) / std::pow( TScalar( 1 ) + A, TScalar( 4 ) ) );
  58     }
  59 };
  60
  61 // -------------------------------------------------------------------------
  62 int main( int argc, char* argv[] )
  63 {
  64   if( argc < 3 )
  65   {
  66     std::cerr
  67       << "Usage: " << argv[ 0 ]
  68       << " input_image neighborhood_order [stop_at_one_front]"
  69       << std::endl;
  70     return( 1 );
  71
  72   } // fi
  73   std::string input_image_fn = argv[ 1 ];
  74   unsigned int neighborhood_order = std::atoi( argv[ 2 ] );
  75   bool stop_at_one_front = false;
  76   if( 3 < argc )
  77     stop_at_one_front = ( std::atoi( argv[ 3 ] ) == 1 );
  78
  79   // Read image
  80   TImageReader::Pointer input_image_reader = TImageReader::New( );
  81   input_image_reader->SetFileName( input_image_fn );
  82   try
  83   {
  84     input_image_reader->Update( );
  85   }
  86   catch( itk::ExceptionObject& err )
  87   {
  88     std::cerr << "Error caught: " << err << std::endl;
  89     return( 1 );
  90
  91   } // yrt
  92   TImage::ConstPointer input_image = input_image_reader->GetOutput( );
  93
  94   TVTKImage::Pointer vtk_image = TVTKImage::New( );
  95   vtk_image->SetInput( input_image );
  96   vtk_image->Update( );
  97
  98   // VTK visualization
  99   vtkSmartPointer< vtkImageActor > actor =
 100     vtkSmartPointer< vtkImageActor >::New( );
 101   actor->SetInputData( vtk_image->GetOutput( ) );
 102
 103   vtkSmartPointer< vtkRenderer > renderer =
 104     vtkSmartPointer< vtkRenderer >::New( );
 105   renderer->SetBackground( 0.1, 0.2, 0.7 );
 106   renderer->AddActor( actor );
 107   vtkSmartPointer< vtkRenderWindow > window =
 108     vtkSmartPointer< vtkRenderWindow >::New( );
 109   window->SetSize( 800, 800 );
 110   window->AddRenderer( renderer );
 111
 112   // VTK interaction
 113   vtkSmartPointer< vtkInteractorStyleImage > imageStyle =
 114     vtkSmartPointer< vtkInteractorStyleImage >::New( );
 115   vtkSmartPointer< vtkRenderWindowInteractor > interactor =
 116     vtkSmartPointer< vtkRenderWindowInteractor >::New( );
 117   interactor->SetInteractorStyle( imageStyle );
 118   window->SetInteractor( interactor );
 119   window->Render( );
 120
 121   // Create the widget and its representation
 122   vtkSmartPointer< vtkPointHandleRepresentation3D > handle =
 123     vtkSmartPointer< vtkPointHandleRepresentation3D >::New( );
 124   handle->GetProperty( )->SetColor( 1, 0, 0 );
 125   vtkSmartPointer< vtkSeedRepresentation > rep =
 126     vtkSmartPointer< vtkSeedRepresentation >::New( );
 127   rep->SetHandleRepresentation( handle );
 128
 129   vtkSmartPointer< vtkSeedWidget > widget =
 130     vtkSmartPointer< vtkSeedWidget >::New( );
 131   widget->SetInteractor( interactor );
 132   widget->SetRepresentation( rep );
 133
 134   // Let some interaction
 135   interactor->Initialize( );
 136   window->Render( );
 137   widget->On( );
 138   interactor->Start( );
 139
 140   // Compute cost map
 141   typedef itk::MinimumMaximumImageCalculator< TImage > TMinMax;
 142   typedef itk::InvertIntensityImageFilter< TImage > TInvert;
 143   typedef itk::DanielssonDistanceMapImageFilter< TImage, TScalarImage > TDanielsson;
 144   typedef itk::UnaryFunctorImageFilter< TScalarImage, TScalarImage, InvertPixelFunctor > TInvertFunctor;
 145
 146   TMinMax::Pointer input_image_min_max = TMinMax::New( );
 147   input_image_min_max->SetImage( input_image );
 148   input_image_min_max->Compute( );
 149
 150   TInvert::Pointer invert = TInvert::New( );
 151   invert->SetInput( input_image );
 152   invert->SetMaximum( input_image_min_max->GetMaximum( ) );
 153
 154   TDanielsson::Pointer danielsson = TDanielsson::New( );
 155   danielsson->SetInput( invert->GetOutput( ) );
 156   danielsson->InputIsBinaryOn( );
 157   danielsson->SquaredDistanceOff( );
 158   danielsson->UseImageSpacingOn( );
 159
 160   TInvertFunctor::Pointer invert_pixels = TInvertFunctor::New( );
 161   invert_pixels->SetInput( danielsson->GetOutput( ) );
 162   invert_pixels->Update( );
 163
 164   // Configure observer
 165   TObserver::Pointer obs = TObserver::New( );
 166   obs->SetImage( invert_pixels->GetOutput( ), window );
 167
 168   // Configure algorithm
 169   TFrontAlgorithm::Pointer algorithm = TFrontAlgorithm::New( );
 170   for( unsigned int s = 0; s < rep->GetNumberOfSeeds( ); s++ )
 171   {
 172     double pos[ 3 ];
 173     rep->GetSeedWorldPosition( s, pos );
 174
 175     TImage::PointType pnt;
 176     pnt[ 0 ] = pos[ 0 ];
 177     pnt[ 1 ] = pos[ 1 ];
 178
 179     TImage::IndexType idx;
 180     if( input_image->TransformPhysicalPointToIndex( pnt, idx ) )
 181     {
 182       algorithm->AddSeed( idx, 0 );
 183       std::cout << " Seed --> " << idx << std::endl;
 184
 185     } // fi
 186
 187   } // rof
 188   algorithm->AddObserver( itk::AnyEvent( ), obs );
 189   algorithm->ThrowEventsOn( );
 190   algorithm->SetInput( invert_pixels->GetOutput( ) );
 191   algorithm->SetNeighborhoodOrder( neighborhood_order );
 192   algorithm->SetStopAtOneFront( stop_at_one_front );
 193   algorithm->Update( );
 194
 195   // One last interaction
 196   window->Render( );
 197   interactor->Start( );
 198
 199   return( 0 );
 200 }
 201
 202 // eof - $RCSfile$