tools/clitkAverageTemporalDimensionGenericFilter.txx

   1 #ifndef clitkAverageTemporalDimensionGenericFilter_txx
   2 #define clitkAverageTemporalDimensionGenericFilter_txx
   3
   4 /* =================================================
   5  * @file   clitkAverageTemporalDimensionGenericFilter.txx
   6  * @author
   7  * @date
   8  *
   9  * @brief
  10  *
  11  ===================================================*/
  12
  13
  14 namespace clitk
  15 {
  16
  17   //-----------------------------------------------------------
  18   // Constructor
  19   //-----------------------------------------------------------
  20   template<class args_info_type>
  21   AverageTemporalDimensionGenericFilter<args_info_type>::AverageTemporalDimensionGenericFilter()
  22   {
  23     m_Verbose=false;
  24     m_InputFileName="";
  25   }
  26
  27
  28   //-----------------------------------------------------------
  29   // Update
  30   //-----------------------------------------------------------
  31   template<class args_info_type>
  32   void AverageTemporalDimensionGenericFilter<args_info_type>::Update()
  33   {
  34     // Read the Dimension and PixelType
  35     int Dimension, Components;
  36     std::string PixelType;
  37     ReadImageDimensionAndPixelType(m_InputFileName, Dimension, PixelType, Components);
  38
  39
  40     // Call UpdateWithDim
  41     if (m_ArgsInfo.input_given>1) Dimension+=1;
  42     if(Dimension==2) UpdateWithDim<2>(PixelType, Components);
  43     else if(Dimension==3)UpdateWithDim<3>(PixelType, Components);
  44     else if (Dimension==4)UpdateWithDim<4>(PixelType, Components);
  45     else
  46       {
  47         std::cout<<"Error, Only for 2, 3 or 4D!!!"<<std::endl ;
  48         return;
  49       }
  50   }
  51
  52   //-------------------------------------------------------------------
  53   // Update with the number of dimensions
  54   //-------------------------------------------------------------------
  55   template<class args_info_type>
  56   template<unsigned int Dimension>
  57   void
  58   AverageTemporalDimensionGenericFilter<args_info_type>::UpdateWithDim(const std::string PixelType, const int Components)
  59   {
  60     if (m_Verbose) std::cout << "Image was detected to be "<<Dimension<<"D with "<<Components<<" component(s) of "<<  PixelType<<"..."<<std::endl;
  61
  62     if (Components==1)
  63       {
  64         if(PixelType == "short"){
  65           if (m_Verbose) std::cout << "Launching filter in "<< Dimension <<"D and signed short..." << std::endl;
  66           UpdateWithDimAndPixelType<Dimension, signed short>();
  67         }
  68         //    else if(PixelType == "unsigned_short"){
  69         //       if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and unsigned_short..." << std::endl;
  70         //       UpdateWithDimAndPixelType<Dimension, unsigned short>();
  71         //     }
  72
  73         else if (PixelType == "unsigned_char"){
  74           if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and unsigned_char..." << std::endl;
  75           UpdateWithDimAndPixelType<Dimension, unsigned char>();
  76         }
  77
  78         //     else if (PixelType == "char"){
  79         //       if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and signed_char..." << std::endl;
  80         //       UpdateWithDimAndPixelType<Dimension, signed char>();
  81         //     }
  82         else {
  83           if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and float..." << std::endl;
  84           UpdateWithDimAndPixelType<Dimension, float>();
  85         }
  86       }
  87
  88     //     else if (Components==2)
  89     //       {
  90     //  if (m_Verbose) std::cout  << "Launching transform in "<< Dimension <<"D and 3D float (DVF)" << std::endl;
  91     //  UpdateWithDimAndPixelType<Dimension, itk::Vector<float, 2> >();
  92     //       }
  93
  94     else if (Components==3)
  95       {
  96         if (m_Verbose) std::cout  << "Launching transform in "<< Dimension <<"D and 3D float (DVF)" << std::endl;
  97         UpdateWithDimAndPixelType<Dimension, itk::Vector<float, 3> >();
  98       }
  99
 100     //     else if (Components==4)
 101     //       {
 102     //  if (m_Verbose) std::cout  << "Launching transform in "<< Dimension <<"D and 3D float (DVF)" << std::endl;
 103     //  UpdateWithDimAndPixelType<Dimension, itk::Vector<float, 3> >();
 104     //       }
 105     else std::cerr<<"Number of components is "<<Components<<", not supported!"<<std::endl;
 106   }
 107
 108
 109   //-------------------------------------------------------------------
 110   // Update with the number of dimensions and the pixeltype
 111   //-------------------------------------------------------------------
 112   template<class args_info_type>
 113   template <unsigned int Dimension, class  PixelType>
 114   void
 115   AverageTemporalDimensionGenericFilter<args_info_type>::UpdateWithDimAndPixelType()
 116   {
 117
 118     // ImageTypes
 119     typedef itk::Image<PixelType, Dimension> InputImageType;
 120     typedef itk::Image<PixelType, Dimension-1> OutputImageType;
 121     typename InputImageType::Pointer input;
 122
 123     if (m_ArgsInfo.input_given ==1 )
 124       {
 125         // Read the input
 126         typedef itk::ImageFileReader<InputImageType> InputReaderType;
 127         typename InputReaderType::Pointer reader = InputReaderType::New();
 128         reader->SetFileName( m_InputFileName);
 129         reader->Update();
 130         input= reader->GetOutput();
 131       }
 132
 133     else
 134       {
 135         // Read and join multiple inputs
 136         if (m_Verbose) std::cout<<m_ArgsInfo.input_given<<" inputs given..."<<std::endl;
 137         std::vector<std::string> filenames;
 138         for(unsigned int i=0; i<m_ArgsInfo.input_given;i++)
 139           {
 140             if (m_Verbose) std::cout<<m_ArgsInfo.input_arg[i]<<std::endl;
 141             filenames.push_back(m_ArgsInfo.input_arg[i]);
 142           }
 143
 144         typedef itk::ImageSeriesReader<InputImageType> ImageReaderType;
 145         typename  ImageReaderType::Pointer reader= ImageReaderType::New();
 146         reader->SetFileNames(filenames);
 147         reader->Update();
 148         input =reader->GetOutput();
 149       }
 150
 151
 152     // Output properties
 153     typename OutputImageType::RegionType region;
 154     typename OutputImageType::RegionType::SizeType size;
 155     typename OutputImageType::IndexType index;
 156     typename OutputImageType::SpacingType spacing;
 157     typename OutputImageType::PointType origin;
 158     typename InputImageType::RegionType region4D=input->GetLargestPossibleRegion();
 159     typename InputImageType::RegionType::SizeType size4D=region4D.GetSize();
 160     typename InputImageType::IndexType index4D=region4D.GetIndex();
 161     typename InputImageType::SpacingType spacing4D=input->GetSpacing();
 162     typename InputImageType::PointType origin4D=input->GetOrigin();
 163
 164     for (unsigned int i=0; i< Dimension-1; i++)
 165       {
 166         size[i]=size4D[i];
 167         index[i]=index4D[i];
 168         spacing[i]=spacing4D[i];
 169         origin[i]=origin4D[i];
 170       }
 171     region.SetSize(size);
 172     region.SetIndex(index);
 173     typename OutputImageType::Pointer output= OutputImageType::New();
 174     output->SetRegions(region);
 175     output->SetSpacing(spacing);
 176     output->SetOrigin(origin);
 177     output->Allocate();
 178
 179
 180     // Region iterators
 181     typedef itk::ImageRegionIterator<InputImageType> IteratorType;
 182     std::vector<IteratorType> iterators(size4D[Dimension-1]);
 183     for (unsigned int i=0; i< size4D[Dimension-1]; i++)
 184       {
 185         typename InputImageType::RegionType regionIt=region4D;
 186         typename InputImageType::RegionType::SizeType sizeIt=regionIt.GetSize();
 187         sizeIt[Dimension-1]=1;
 188         regionIt.SetSize(sizeIt);
 189         typename InputImageType::IndexType indexIt=regionIt.GetIndex();
 190         indexIt[Dimension-1]=i;
 191         regionIt.SetIndex(indexIt);
 192         iterators[i]=IteratorType(input, regionIt);
 193        }
 194
 195     typedef itk::ImageRegionIterator<OutputImageType> OutputIteratorType;
 196     OutputIteratorType avIt(output, output->GetLargestPossibleRegion());
 197
 198     // Average
 199     PixelType vector;
 200     PixelType zeroVector=itk::NumericTraits<PixelType>::Zero;
 201     //zeroVector.Fill(0.0);
 202     while (!(iterators[0]).IsAtEnd())
 203       {
 204         vector=zeroVector;
 205         for (unsigned int i=0; i<size4D[Dimension-1]; i++)
 206           {
 207             vector+=iterators[i].Get();
 208             ++(iterators[i]);
 209           }
 210         vector/=size4D[Dimension-1];
 211         avIt.Set(vector);
 212         ++avIt;
 213       }
 214
 215     // Output
 216     typedef itk::ImageFileWriter<OutputImageType> WriterType;
 217     typename WriterType::Pointer writer = WriterType::New();
 218     writer->SetFileName(m_ArgsInfo.output_arg);
 219     writer->SetInput(output);
 220     writer->Update();
 221
 222   }
 223
 224
 225 }//end clitk
 226
 227 #endif //#define clitkAverageTemporalDimensionGenericFilter_txx