appli/CTBronchi/RandomWalker.cxx

   1 #include <chrono>
   2 #include <itkImage.h>
   3 #include <itkImageFileReader.h>
   4 #include <itkImageFileWriter.h>
   5
   6 #include <fpa/Image/RandomWalker.h>
   7 #include <fpa/Image/Functors/Dijkstra/Gaussian.h>
   8 #include <ivq/ITK/ExtractLabelFunction.h>
   9 #include <ivq/ITK/ImageUnaryFunctionFilter.h>
  10
  11 // -------------------------------------------------------------------------
  12 const unsigned int Dim = 3;
  13 typedef short          TPixel;
  14 typedef unsigned short TLabel;
  15 typedef unsigned char  TBinary;
  16 typedef double         TScalar;
  17
  18 typedef itk::Image< TPixel, Dim > TInputImage;
  19 typedef itk::Image< TLabel, Dim > TLabelImage;
  20 typedef itk::Image< TBinary, Dim > TBinaryImage;
  21
  22 typedef fpa::Image::RandomWalker< TInputImage, TLabelImage, TScalar > TFilter;
  23 typedef fpa::Image::Functors::Dijkstra::Gaussian< TInputImage, TScalar > TWeight;
  24
  25 typedef ivq::ITK::ExtractLabelFunction< TLabel, TBinary > TLabelFunction;
  26 typedef ivq::ITK::ImageUnaryFunctionFilter< TLabelImage, TBinaryImage > TLabelExtractor;
  27
  28 // -------------------------------------------------------------------------
  29 int main( int argc, char* argv[] )
  30 {
  31   // Get arguments
  32   if( argc < 7 )
  33   {
  34     std::cerr
  35       << "Usage: " << argv[ 0 ] << std::endl
  36       << "   input_image label_image output_image" << std::endl
  37       << "   label alpha(0) beta(100)"
  38       << std::endl;
  39     return( 1 );
  40
  41   } // fi
  42   std::string input_image_filename = argv[ 1 ];
  43   std::string label_image_filename = argv[ 2 ];
  44   std::string output_image_filename = argv[ 3 ];
  45   TLabel label = TLabel( std::atoi( argv[ 4 ] ) );
  46   double alpha = std::atof( argv[ 5 ] );
  47   double beta = std::atof( argv[ 6 ] );
  48
  49   // Read images
  50   itk::ImageFileReader< TInputImage >::Pointer input_image_reader =
  51     itk::ImageFileReader< TInputImage >::New( );
  52   input_image_reader->SetFileName( input_image_filename );
  53
  54   itk::ImageFileReader< TLabelImage >::Pointer label_image_reader =
  55     itk::ImageFileReader< TLabelImage >::New( );
  56   label_image_reader->SetFileName( label_image_filename );
  57
  58   // Prepare weight
  59   TWeight::Pointer weight = TWeight::New( );
  60   weight->SetAlpha( alpha );
  61   weight->SetBeta( beta );
  62
  63   // Prepare filter
  64   TFilter::Pointer filter = TFilter::New( );
  65   filter->SetInput( input_image_reader->GetOutput( ) );
  66   filter->SetLabels( label_image_reader->GetOutput( ) );
  67   filter->SetWeightFunction( weight );
  68
  69   // Label extraction
  70   TLabelFunction::Pointer label_function = TLabelFunction::New( );
  71   label_function->SetLabel( label );
  72   label_function->SetInsideValue( 255 );
  73   label_function->SetOutsideValue( 0 );
  74
  75   TLabelExtractor::Pointer label_extractor = TLabelExtractor::New( );
  76   label_extractor->SetInput( filter->GetMarks( ) );
  77   label_extractor->SetFunction( label_function );
  78
  79
  80   // Show some information
  81   std::cout << "----------------------------------------------" << std::endl;
  82   std::cout << "Image: " << input_image_filename << std::endl;
  83
  84   // Execute pipeline
  85   std::chrono::time_point< std::chrono::high_resolution_clock > ts, te;
  86   std::chrono::duration< double > tr;
  87   try
  88   {
  89     ts = std::chrono::high_resolution_clock::now( );
  90     input_image_reader->Update( );
  91     te = std::chrono::high_resolution_clock::now( );
  92     tr = te - ts;
  93     std::cout << "Raw read time: " << tr.count( ) << " s" << std::endl;
  94
  95     ts = std::chrono::high_resolution_clock::now( );
  96     label_image_reader->Update( );
  97     te = std::chrono::high_resolution_clock::now( );
  98     tr = te - ts;
  99     std::cout << "Label read time: " << tr.count( ) << " s" << std::endl;
 100
 101     ts = std::chrono::high_resolution_clock::now( );
 102     filter->Update( );
 103     te = std::chrono::high_resolution_clock::now( );
 104     tr = te - ts;
 105     std::cout
 106       << "Labelling time: " << tr.count( ) << " s" << std::endl;
 107
 108     ts = std::chrono::high_resolution_clock::now( );
 109     label_extractor->Update( );
 110     te = std::chrono::high_resolution_clock::now( );
 111     tr = te - ts;
 112     std::cout
 113       << "Extract label time: " << tr.count( ) << " s" << std::endl;
 114   }
 115   catch( std::exception& err )
 116   {
 117     std::cerr << "Error caught: " << err.what( ) << std::endl;
 118     return( 1 );
 119
 120   } // yrt
 121
 122   // Save output image
 123   itk::ImageFileWriter< TBinaryImage >::Pointer output_image_writer =
 124     itk::ImageFileWriter< TBinaryImage >::New( );
 125   output_image_writer->SetInput( label_extractor->GetOutput( ) );
 126   output_image_writer->SetFileName( output_image_filename );
 127   try
 128   {
 129     ts = std::chrono::high_resolution_clock::now( );
 130     output_image_writer->Update( );
 131     te = std::chrono::high_resolution_clock::now( );
 132     tr = te - ts;
 133     std::cout << "Write time: " << tr.count( ) << " s" << std::endl;
 134   }
 135   catch( std::exception& err )
 136   {
 137     std::cerr << "Error caught: " << err.what( ) << std::endl;
 138     return( 1 );
 139
 140   } // yrt
 141   std::cout << "----------------------------------------------" << std::endl;
 142
 143   return( 0 );
 144 }
 145
 146 // eof - $RCSfile$