[Insight-users] ITK example that sums multiple metaheader images

[Stéphane CALANDE]

Hi.

I just send this for anyone that could be interrested in the
itkNaryMeanImageFilter created for my problem.

There was a little problem of types.
With the code of "itkNaryMeanImageFilter" that was given in the former mail,
if you wanted to calculate the mean of a lot of images, there was an
overflow of the variable "AccumulatorType mean = NumericTraits< TOutput
>::Zero;" in the "for" :

for( unsigned int i=0; i< input.size(); i++ )
     {
     mean += static_cast< TOutput >(input[i])
     }

Even I divied "static_cast< TOutput >(input[i])" by the number of images in
the "for", I had till problems (some parts of the image were black)

Thats why I replace that part of the code by that one :


//AccumulatorType mean = NumericTraits< TOutput >::Zero;
long double total = 0;
   for( unsigned int i=0; i< input.size(); i++ )
     {
     //mean += (static_cast< TOutput >(input[i])) / input.size() ;
     total += static_cast< TOutput >(input[i]) ;
     }
   //return static_cast<TOutput>( mean );
   total = total / input.size();
   return static_cast<TOutput>( total );



I give you here the final code of itkNaryMeanImageFilter that works + my
program that take multiple mhd images and that give as output the mean of
that images (Sorry the "cout" are in French ;-) )

I hope it will be useful for anyone interrested...



Stéphane.


/*=========================================================================

 Program:   Insight Segmentation & Registration Toolkit
 Module:    $RCSfile: itkNaryMeanImageFilter.h,v $
 Language:  C++
 Date:      $Date$
 Version:   $Revision$

 Copyright (c) Insight Software Consortium. All rights reserved.
 See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.

    This software is distributed WITHOUT ANY WARRANTY; without even
    the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
    PURPOSE.  See the above copyright notices for more information.

=========================================================================*/
#ifndef __itkNaryMeanImageFilter_h
#define __itkNaryMeanImageFilter_h

#include "itkNaryFunctorImageFilter.h"
#include "itkNumericTraits.h"

namespace itk
{

/** \class NaryMeanImageFilter
 * \brief Implements an operator for pixel-wise mean of N images.
 *
 * This class is parametrized over the types of the two
 * input images and the type of the output image.
 * Numeric conversions (castings) are done by the C++ defaults.
 *
 * The pixel type of the input 1 image must have a valid defintion of
 * the operator+ with a pixel type of the image 2. This condition is
 * required because internally this filter will perform the operation
 *
 *        sum of input pixels / number of inputs
 *
 * Additionally the type resulting from the mean, will be cast to
 * the pixel type of the output image.
 *
 * The total operation over one pixel will be
 *
 *  output_pixel = static_cast<OutputPixelType>( sum of input pixels /
num inputs )
 *
 * \warning No numeric overflow checking is performed in this filter.
 *
 * \ingroup IntensityImageFilters  Multithreaded
 */

namespace Functor {

template< class TInput, class TOutput >
class Mean
{
public:
 typedef typename NumericTraits< TInput >::AccumulateType AccumulatorType;
 Mean() {}
 ~Mean() {}
 inline TOutput operator()( const std::vector< TInput > & input )
 {
   long double total = 0;
   for( unsigned int i=0; i< input.size(); i++ )
     {
     total += static_cast< TOutput >(input[i]) ;
     }
   total = total / input.size();
   return static_cast<TOutput>( total );
 }
 bool operator== (const Mean&) const
 {
   return true;
 }
 bool operator!= (const Mean&) const
 {
   return false;
 }
};
}
template <class TInputImage, class TOutputImage>
class ITK_EXPORT NaryMeanImageFilter :
   public
NaryFunctorImageFilter<TInputImage,TOutputImage,
                      Functor::Mean<typename TInputImage::PixelType,
typename TInputImage::PixelType > >
{
public:
 /** Standard class typedefs. */
 typedef NaryMeanImageFilter  Self;
 typedef NaryFunctorImageFilter<TInputImage,TOutputImage,
                                Functor::Mean<typename
TInputImage::PixelType,
                                typename TInputImage::PixelType  > >
Superclass;
 typedef SmartPointer<Self>   Pointer;
 typedef SmartPointer<const Self>  ConstPointer;

 /** Method for creation through the object factory. */
 itkNewMacro(Self);

 /** Runtime information support. */
 itkTypeMacro(NaryMeanImageFilter,
              NaryFunctorImageFilter);

#ifdef ITK_USE_CONCEPT_CHECKING
 /** Begin concept checking */
 itkConceptMacro(InputConvertibleToOutputCheck,
   (Concept::Convertible<typename TInputImage::PixelType,
                         typename TOutputImage::PixelType>));
 itkConceptMacro(InputHasZeroCheck,
   (Concept::HasZero<typename TInputImage::PixelType>));
 /** End concept checking */
#endif

protected:
 NaryMeanImageFilter() {}
 virtual ~NaryMeanImageFilter() {}

private:
 NaryMeanImageFilter(const Self&); //purposely not implemented
 void operator=(const Self&); //purposely not implemented

};

} // end namespace itk


#endif




**********************************************************************************************************************************************
**********************************************************************************************************************************************
**********************************************************************************************************************************************



 #include "itkImage.h"
 #include "itkImageFileReader.h"
 #include "itkImageFileWriter.h"
 #include "itkImageIOBase.h"
 #include "itkNaryMeanImageFilter.h"

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

 {

 if( argc < 4 )

 {

 std::cerr << "Usage: " << std::endl;

 std::cerr << argv[0] << " Image1.mhd Image2.mhd [ImageX.mhd]*
NomImageOutput.mhd" << std::endl;

 return EXIT_FAILURE;

 }

 typedef   signed short PixelType; // être sûr que c'est le bon type !!

 typedef itk::Image< PixelType, 3 > ImageType;

 typedef itk::ImageFileReader< ImageType > ReaderType;

 typedef itk::ImageFileWriter< ImageType > WriterType;

 typedef itk::NaryMeanImageFilter< ImageType,
                                  ImageType > MeanFilterType;

 int nbImages = argc - 2;

 ReaderType::Pointer reader = ReaderType::New();

 MeanFilterType::Pointer addition = MeanFilterType::New();

 std::cout << std::endl;

 for (int i = 0 ; i < nbImages ; i++)
{
    reader->SetFileName( argv[i+1] );
    reader->Update();
    ImageType::Pointer input = reader->GetOutput();
    input->DisconnectPipeline();
    addition->SetInput( i , input );
    std::cout << "Image n°" << i+1 << " [" << argv[i+1] << "] ajoutée" <<
std::endl << std::endl;
}

 addition-> Update();

 WriterType::Pointer writer = WriterType::New();

 writer->SetFileName( argv[nbImages+1] );

 writer->SetInput(addition->GetOutput());

 std::cout << "Ecriture du fichier..." << std::endl << std::endl;

 writer->Update();

 std::cout << "Fichier '" << argv[nbImages+1] << "' créé" << std::endl;

 return EXIT_SUCCESS;

 }
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20081024/56cac135/attachment-0001.htm>