tests/image/RandomWalker/CompareRandomWalkers.cxx

   1 // =========================================================================
   2 // @author Leonardo Florez Valencia
   3 // @email florez-l@javeriana.edu.co
   4 // =========================================================================
   5
   6 #include <chrono>
   7 #include <iostream>
   8 #include <string>
   9 #include <sstream>
  10 #include <itkImage.h>
  11 #include <itkImageFileReader.h>
  12 #include <itkImageFileWriter.h>
  13
  14 #include <fpa/Common/OriginalRandomWalker.h>
  15 #include <fpa/Filters/Image/RandomWalker.h>
  16 #include <fpa/Functors/Dijkstra/Image/Gaussian.h>
  17
  18 // -------------------------------------------------------------------------
  19 const unsigned int Dim = 3;
  20 typedef short         TPixel;
  21 typedef unsigned char TLabel;
  22 typedef double        TScalar;
  23
  24 typedef itk::Image< TPixel, Dim > TImage;
  25 typedef itk::Image< TLabel, Dim > TLabels;
  26 typedef itk::Image< TScalar, Dim > TScalarImage;
  27
  28 // -------------------------------------------------------------------------
  29 double MeasureTime( itk::ProcessObject* f )
  30 {
  31   std::chrono::time_point< std::chrono::high_resolution_clock > s, e;
  32   std::chrono::duration< double > t;
  33   s = std::chrono::high_resolution_clock::now( );
  34   f->Update( );
  35   e = std::chrono::high_resolution_clock::now( );
  36   t = e - s;
  37   return( t.count( ) );
  38 }
  39
  40 // -------------------------------------------------------------------------
  41 template< class _TImagePtr >
  42 void ReadImage( _TImagePtr& image, const std::string& fname )
  43 {
  44   typedef itk::ImageFileReader< typename _TImagePtr::ObjectType > _TReader;
  45   typename _TReader::Pointer r = _TReader::New( );
  46   r->SetFileName( fname );
  47   double t = MeasureTime( r );
  48   image = r->GetOutput( );
  49   image->DisconnectPipeline( );
  50   std::cout << "Image \"" << fname << "\" read in " << t << "s" << std::endl;
  51 }
  52
  53 // -------------------------------------------------------------------------
  54 template< class _TImagePtr >
  55 void WriteImage( const _TImagePtr& image, const std::string& fname )
  56 {
  57   typedef itk::ImageFileWriter< typename _TImagePtr::ObjectType > _TWriter;
  58   typename _TWriter::Pointer w = _TWriter::New( );
  59   w->SetInput( image );
  60   w->SetFileName( fname );
  61   double t = MeasureTime( w );
  62   std::cout << "Image \"" << fname << "\" written in " << t << "s" << std::endl;
  63 }
  64
  65 // -------------------------------------------------------------------------
  66 template< class _TImagePtr, class _TLabelsPtr, class _TScalarImagePtr >
  67 void Original(
  68   _TLabelsPtr& output_labels,
  69   _TScalarImagePtr& output_probabilities,
  70   const _TImagePtr& input,
  71   const _TLabelsPtr& labels,
  72   const typename _TScalarImagePtr::ObjectType::PixelType& beta,
  73   const typename _TScalarImagePtr::ObjectType::PixelType& epsilon
  74   )
  75 {
  76   typedef typename _TImagePtr::ObjectType       _TImage;
  77   typedef typename _TLabelsPtr::ObjectType      _TLabels;
  78   typedef typename _TScalarImagePtr::ObjectType _TScalarImage;
  79   typedef typename _TScalarImage::PixelType     _TScalar;
  80   typedef
  81     fpa::Functors::Dijkstra::Image::Gaussian< _TImage, _TScalar >
  82     _TFunction;
  83   typedef
  84     fpa::Common::OriginalRandomWalker< _TImage, _TLabels, _TScalar >
  85     _TFilter;
  86
  87   // Random walker function
  88   typename _TFunction::Pointer edge = _TFunction::New( );
  89   edge->SetBeta( beta );
  90   edge->SetEpsilon( epsilon );
  91   edge->TreatAsWeightOff( );
  92
  93   // Random walker
  94   typename _TFilter::Pointer filter = _TFilter::New( );
  95   filter->SetInput( input );
  96   filter->SetInputLabels( labels );
  97   filter->SetEdgeFunction( edge );
  98   filter->DebugOn( );
  99   double t = MeasureTime( filter );
 100   output_labels = filter->GetOutput( );
 101   output_probabilities = filter->GetOutputProbabilities( );
 102   output_labels->DisconnectPipeline( );
 103   output_probabilities->DisconnectPipeline( );
 104   std::cout << "Original RW executed in " << t << "s" << std::endl;
 105 }
 106
 107 // -------------------------------------------------------------------------
 108 template< class _TImagePtr, class _TLabelsPtr, class _TScalarImagePtr >
 109 void FrontPropagation(
 110   _TLabelsPtr& output_labels,
 111   _TScalarImagePtr& output_probabilities,
 112   const _TImagePtr& input,
 113   const _TLabelsPtr& labels,
 114   const typename _TScalarImagePtr::ObjectType::PixelType& beta,
 115   const typename _TScalarImagePtr::ObjectType::PixelType& epsilon
 116   )
 117 {
 118   typedef typename _TImagePtr::ObjectType       _TImage;
 119   typedef typename _TLabelsPtr::ObjectType      _TLabels;
 120   typedef typename _TScalarImagePtr::ObjectType _TScalarImage;
 121   typedef typename _TScalarImage::PixelType     _TScalar;
 122   typedef
 123     fpa::Functors::Dijkstra::Image::Gaussian< _TImage, _TScalar >
 124     _TFunction;
 125   typedef
 126     fpa::Filters::Image::RandomWalker< _TImage, _TLabels, _TScalar >
 127     _TFilter;
 128
 129   // Random walker function
 130   typename _TFunction::Pointer edge = _TFunction::New( );
 131   edge->SetBeta( beta );
 132   edge->SetEpsilon( epsilon );
 133   edge->TreatAsWeightOn( );
 134
 135   // Random walker
 136   typename _TFilter::Pointer filter = _TFilter::New( );
 137   filter->SetInputImage( input );
 138   filter->SetInputLabels( labels );
 139   filter->SetWeightFunction( edge );
 140   double t = MeasureTime( filter );
 141   output_labels = filter->GetOutputLabels( );
 142   output_probabilities = filter->GetOutputCosts( );
 143   output_labels->DisconnectPipeline( );
 144   output_probabilities->DisconnectPipeline( );
 145   std::cout << "Front propagation RW executed in " << t << "s" << std::endl;
 146 }
 147
 148 // -------------------------------------------------------------------------
 149 int main( int argc, char* argv[] )
 150 {
 151   // Get input arguments
 152   if( argc < 5 )
 153   {
 154     std::cerr
 155       << "Usage: " << argv[ 0 ]
 156       << " input_image input_labels beta epsilon"
 157       << std::endl;
 158     return( 1 );
 159
 160   } // fi
 161   std::string input_image_filename = argv[ 1 ];
 162   std::string input_labels_filename = argv[ 2 ];
 163   std::istringstream values( std::string( argv[ 3 ] ) + " " + argv[ 4 ] );
 164   TScalar beta, epsilon;
 165   values >> beta >> epsilon;
 166
 167   // Read inputs
 168   TImage::Pointer input_image;
 169   TLabels::Pointer input_labels;
 170   ReadImage( input_image, input_image_filename );
 171   ReadImage( input_labels, input_labels_filename );
 172
 173   // Execute algorithms
 174   TLabels::Pointer original_labels, fp_labels;
 175   TScalarImage::Pointer original_probabilities, fp_probabilities;
 176   Original(
 177     original_labels, original_probabilities,
 178     input_image, input_labels,
 179     beta, epsilon
 180     );
 181   FrontPropagation(
 182     fp_labels, fp_probabilities,
 183     input_image, input_labels,
 184     beta, epsilon
 185     );
 186
 187   WriteImage( original_labels, "orig.mhd" );
 188   WriteImage( fp_labels, "fp.mhd" );
 189
 190   return( 0 );
 191 }
 192
 193 // eof - $RCSfile$