appli/examples/example_ImageGaussianModelEstimator.cxx

   1 #include <cstdlib>
   2 #include <iostream>
   3 #include <string>
   4
   5 #include <itkImage.h>
   6 #include <itkImageFileReader.h>
   7 #include <itkImageRegionConstIterator.h>
   8 #include <itkMatrix.h>
   9 #include <itkRGBPixel.h>
  10 #include <itkVector.h>
  11
  12 #include <cpPlugins/Extensions/Algorithms/IterativeGaussianModelEstimator.h>
  13
  14 // -------------------------------------------------------------------------
  15 const unsigned int Dim = 2;
  16 typedef double                    TScalar;
  17 typedef unsigned char             TChannel;
  18 typedef itk::RGBPixel< TChannel > TPixel;
  19 typedef itk::Image< TPixel, Dim > TImage;
  20
  21 // -------------------------------------------------------------------------
  22 int main( int argc, char* argv[] )
  23 {
  24   if( argc < 2 )
  25   {
  26     std::cerr
  27       << "Usage: " << argv[ 0 ]
  28       << " input_image" << std::endl;
  29     return( 1 );
  30
  31   } // fi
  32   std::string input_image_file = argv[ 1 ];
  33
  34   // Read image
  35   itk::ImageFileReader< TImage >::Pointer reader =
  36     itk::ImageFileReader< TImage >::New( );
  37   reader->SetFileName( input_image_file );
  38   try
  39   {
  40     reader->Update( );
  41   }
  42   catch( itk::ExceptionObject& err )
  43   {
  44     std::cerr << "Error caught: " << err << std::endl;
  45     return( 1 );
  46
  47   } // yrt
  48   TImage::ConstPointer input_image = reader->GetOutput( );
  49
  50   // Compute gaussian model
  51   typedef cpPlugins::Extensions::Algorithms::IterativeGaussianModelEstimator< TScalar, 3 > TEstimator;
  52
  53   TEstimator::Pointer estimator = TEstimator::New( );
  54   estimator->Clear( );
  55
  56   itk::ImageRegionConstIterator< TImage > it(
  57     input_image,
  58     input_image->GetRequestedRegion( )
  59     );
  60   for( it.GoToBegin( ); !it.IsAtEnd( ); ++it )
  61   {
  62     estimator->AddSample(
  63       TScalar( it.Get( ).GetRed( ) ),
  64       TScalar( it.Get( ).GetGreen( ) ),
  65       TScalar( it.Get( ).GetBlue( ) )
  66       );
  67
  68   } // rof
  69
  70   itk::Vector< TScalar, 3 > mean;
  71   itk::Matrix< TScalar, 3, 3 > cova;
  72
  73   estimator->GetModel( mean, cova );
  74   std::cout
  75     << "Covariance = " << std::endl << cova << std::endl;
  76   std::cout
  77     << "Inverse covariance = "
  78     << std::endl << cova.GetInverse( ) << std::endl;
  79   std::cout
  80     << "Mean = " << std::endl << mean << std::endl;
  81
  82   return( 0 );
  83 }
  84
  85 // eof - $RCSfile$