appli/experiments/Skeleton_00.cxx

   1 #include <itkImage.h>
   2 #include <itkMinimumMaximumImageCalculator.h>
   3 #include <itkSignedMaurerDistanceMapImageFilter.h>
   4 #include <itkImageFileReader.h>
   5 #include <itkImageFileWriter.h>
   6 #include <fpa/Image/SkeletonFilter.h>
   7 #include <cpExtensions/Algorithms/SkeletonToImageFilter.h>
   8
   9 typedef double TScalar;
  10 typedef itk::Image< TScalar, 3 > TImage;
  11 typedef itk::MinimumMaximumImageCalculator< TImage > TMinMax;
  12 typedef itk::SignedMaurerDistanceMapImageFilter< TImage, TImage > TDMap;
  13 typedef fpa::Image::SkeletonFilter< TImage > TSkeletonFilter;
  14 typedef TSkeletonFilter::TSkeleton TSkeleton;
  15 typedef cpExtensions::Algorithms::SkeletonToImageFilter< TSkeleton, TImage > TSkeletonToImage;
  16 typedef itk::ImageFileReader< TImage > TReader;
  17 typedef itk::ImageFileWriter< TImage > TWriter;
  18
  19 int main( int argc, char* argv[] )
  20 {
  21   if( argc < 4 )
  22   {
  23     std::cerr
  24       << "Usage: " << argv[ 0 ]
  25       << " image prefix number_of_samples [seed_x seed_y seed_z]"
  26       << std::endl;
  27     return( 1 );
  28
  29   } // fi
  30   std::string image_name = argv[ 1 ];
  31
  32   TReader::Pointer reader = TReader::New( );
  33   reader->SetFileName( image_name );
  34   reader->Update( );
  35
  36   TDMap::Pointer dmap = TDMap::New( );
  37   dmap->SetInput( reader->GetOutput( ) );
  38   dmap->InsideIsPositiveOn( );
  39   dmap->Update( );
  40   const TImage* image = dmap->GetOutput( );
  41
  42   std::string prefix = argv[ 2 ];
  43   std::istringstream str( argv[ 3 ] );
  44   unsigned int number_of_samples;
  45   str >> number_of_samples;
  46   TImage::IndexType seed;
  47   if( argc < 7 )
  48   {
  49     TMinMax::Pointer minmax = TMinMax::New( );
  50     minmax->SetImage( image );
  51     minmax->Compute( );
  52     seed = minmax->GetIndexOfMaximum( );
  53   }
  54   else
  55   {
  56     std::istringstream str_x( argv[ 4 ] );
  57     std::istringstream str_y( argv[ 5 ] );
  58     std::istringstream str_z( argv[ 6 ] );
  59     TImage::PointType pnt;
  60     str_x >> pnt[ 0 ];
  61     str_y >> pnt[ 1 ];
  62     str_z >> pnt[ 2 ];
  63     image->TransformPhysicalPointToIndex( pnt, seed );
  64
  65   } // fi
  66
  67   for( unsigned int i = 0; i < number_of_samples; ++i )
  68   {
  69     std::cout << i << std::endl;
  70     TSkeletonFilter::Pointer skeleton = TSkeletonFilter::New( );
  71     skeleton->SetInput( image );
  72     skeleton->AddSeed( seed, 0 );
  73     skeleton->Update( );
  74
  75     TSkeletonToImage::Pointer sk2im = TSkeletonToImage::New( );
  76     sk2im->SetTemplateImage( image );
  77     sk2im->SetSkeleton( skeleton->GetSkeleton( ) );
  78     sk2im->SetOutsideValue( 0 );
  79     sk2im->SetInsideValue( 1 );
  80     sk2im->Update( );
  81
  82     TWriter::Pointer writer = TWriter::New( );
  83     writer->SetInput( sk2im->GetOutput( ) );
  84     writer->SetFileName( "skeleton.mhd" );
  85     writer->Update( );
  86
  87   } // rof
  88   return( 0 );
  89 }
  90
  91 // eof - $RCSfile$