[Insight-users] Problems using SampleMeanShiftBlurringFilter...

[Antonin Perrot-Audet]

Hello All, I am trying to use the Mean Shift Filter. I used the code 
from insightclopedia, but I also need the spacial information.
The final goal would be to use the clustering filter, but first, this is 
important to have the Blurring Filter working and giving good results.
Right now, I have very poor results, and I don't know where it comes 
from : It seems like the mean shift filter is just rescaling the 
intensities of my image !
My idea was that it would identify the blobs as modes ?

I attached the input image I am working on (blobs.png).


int main( int argc, char* argv[] )
{

   //**********TYPEDEF**********
   // Dimension of the image
   const unsigned Dimension = 2;
   // The type of manipulated pixels :
   typedef unsigned char PixelType;

   typedef itk::Vector< unsigned char, 1 >  ProcessingVectorPixelType;


   // The type of manipulated images :
   typedef itk::Image< PixelType, Dimension > ImageType;
   typedef itk::Image< ProcessingVectorPixelType, Dimension > 
VectorImageType;


   //**********ARGUMENTS READING*********
   // We translate command line arguments to be interpreted by our example

   // We do a very basic check : does the user specify enough parameters :
   //  3 inputs ( program path, input_filename, and output_filename )
   if ( argc != 3 )
     {
     std::cerr << "Usage: " << std::endl;
     std::cerr << argv[0] << "(.exe) takes 2 arguments" <<std::endl;
     std::cerr << "1-Input Filename :" <<std::endl;
     std::cerr << "2-Output Filename :" <<std::endl;
     return EXIT_FAILURE;
     }

   // We convert input arguments
   // Input filename
   const char* ArgInputFilename = argv[1]; //from string to string
   const char* ArgOutputFilename = argv[2]; //from string to string




   //**********INPUT IMAGE READING**********

   // A file reader (input) :
   typedef itk::ImageFileReader< VectorImageType > ReaderType;

   ReaderType::Pointer ImageReader = ReaderType::New();  //Create the reader

   ImageReader->SetFileName( ArgInputFilename ); // Set the reader input 
filename

  // try to read the input file
     ImageReader->Update();


   //**********ADAPTING THE INPUT FOR STATISTICAL OPERATIONS**********

   // We store the read image in InputImage :
   VectorImageType::Pointer InputImage = ImageReader->GetOutput();

   // List of Samples Adaptor creation :
   typedef itk::Statistics
              ::JointDomainImageToListAdaptor< VectorImageType > 
ListSampleType;

   //(ScalarImageToListSampleAdaptor in the new statistical framework)
   ListSampleType::Pointer ListSample = ListSampleType::New();

   ListSample->SetImage( InputImage );



     // **********TUNNING THE MEAN SHIFT : SCALE PARAMETERS**********

     ListSampleType::NormalizationFactorsType NormalizationFactors;
     NormalizationFactors [0] = 6; // scale factor along x
     NormalizationFactors [1] = 6; // scale factor along y
     NormalizationFactors [2] = 3; // scale factor along intensities

   ListSample->SetNormalizationFactors(NormalizationFactors);


   // *********ADAPTING THE SAMPLE THE MEAN SHIFT**********
   // Short image samples in a kdTree
   typedef itk::Statistics
              ::KdTreeGenerator< ListSampleType > KdTreeGeneratorType;
   KdTreeGeneratorType::Pointer KdTreeGenerator = 
KdTreeGeneratorType::New() ;
   KdTreeGenerator->SetSample( ListSample ) ;
   KdTreeGenerator->SetBucketSize( 20 ) ;
   KdTreeGenerator->Update() ;

   typedef KdTreeGeneratorType::KdTreeType KdTreeType ;
   KdTreeType::Pointer KdTree = KdTreeGenerator->GetOutput() ;

   typedef itk::Statistics::HypersphereKernelMeanShiftModeSeeker< 
KdTreeType >
              ModeSeekerType;
   ModeSeekerType::Pointer ModeSeeker = ModeSeekerType::New();
   ModeSeeker->SetInputSample( KdTree );
   ModeSeeker->SetSearchRadius(1);


   // Cache for storing temporary results to speedup the modeseeking
   typedef itk::Statistics
              ::MeanShiftModeCacheMethod< 
KdTreeType::MeasurementVectorType >
              CacheMethodType ;
   CacheMethodType::Pointer CacheMethod = CacheMethodType::New();
   CacheMethod->SetMaximumEntries(255);
   CacheMethod->SetMaximumConsecutiveFailures(100);
   CacheMethod->SetHitRatioThreshold( 0.5 );
   ModeSeeker->SetCacheMethod( CacheMethod.GetPointer() );


   //**********MEAN SHIFT**********

   typedef itk::Statistics::SampleMeanShiftBlurringFilter< KdTreeType > 
MSFilterType;
   MSFilterType::Pointer MSFilter = MSFilterType::New();
   MSFilter->SetInputSample( KdTree ) ;
   MSFilter->SetMeanShiftModeSeeker( ModeSeeker );
   MSFilter->Update() ;

   std::cout << "Cache statistics: " << std::endl;
   CacheMethod->Print(std::cout) ;

     //**********MEAN SHIFT OUTPUT TO IMAGE CONVERSION**********

   typedef ImageType OutputImageType ;
   OutputImageType::Pointer OutputImage = OutputImageType::New();
   // The output image is the same size as the input image :
   OutputImage->SetRegions( InputImage->GetLargestPossibleRegion() );
   OutputImage->Allocate() ;

   typedef itk::ImageRegionIterator< OutputImageType > ImageIteratorType;
   ImageIteratorType io_img_iter( OutputImage,
                              OutputImage->GetLargestPossibleRegion() );
   io_img_iter.GoToBegin() ;

   // define iterators over the output of the filter
   // output of the MeanShift Filter
   MSFilterType::OutputType::Pointer MSOutput = MSFilter->GetOutput();

   std::cout << "Mesurement Vector Size (MeanShift Output)"
<< MSOutput->GetMeasurementVectorSize() << std::endl;


   MSFilterType::OutputType::Iterator MSfo_iter = MSOutput->Begin();
   MSFilterType::OutputType::Iterator MSfo_end = MSOutput->End();

   while ( MSfo_iter != MSfo_end )//create output image from measurement 
vector's
                                  //third coordinate : intensity
     {
     io_img_iter.Set( (PixelType) MSfo_iter.GetMeasurementVector()[2]);
     ++MSfo_iter ;
     ++io_img_iter ;
     }


   // *********WRITING OUTPUT IMAGE**********

   // Writer Creation :
   typedef itk::ImageFileWriter< ImageType > WriterType;
   WriterType::Pointer ImageWriter = WriterType::New();

   // Writer initialization (input and output)
   ImageWriter->SetFileName( ArgOutputFilename );  // 
Threshold_filter->writer
   ImageWriter->SetInput( OutputImage);  // writer->Output_filename

   ImageWriter->Update();  // This update will update the whole pipeline



   return EXIT_SUCCESS;
}




Thanks a lot in advance !

-- 
Antonin

-------------- next part --------------
A non-text attachment was scrubbed...
Name: blobs.png
Type: image/png
Size: 29830 bytes
Desc: not available
URL: <http://www.itk.org/pipermail/insight-users/attachments/20100421/db52320d/attachment-0001.png>