tools/clitkBSplineCoefficientsToValuesGenericFilter.txx

   1 /*=========================================================================
   2   Program:   vv                     http://www.creatis.insa-lyon.fr/rio/vv
   3
   4   Authors belong to:
   5   - University of LYON              http://www.universite-lyon.fr/
   6   - Léon Bérard cancer center       http://oncora1.lyon.fnclcc.fr
   7   - CREATIS CNRS laboratory         http://www.creatis.insa-lyon.fr
   8
   9   This software is distributed WITHOUT ANY WARRANTY; without even
  10   the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  11   PURPOSE.  See the copyright notices for more information.
  12
  13   It is distributed under dual licence
  14
  15   - BSD        See included LICENSE.txt file
  16   - CeCILL-B   http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
  17 ======================================================================-====*/
  18 #ifndef clitkBSplineCoefficientsToValuesGenericFilter_txx
  19 #define clitkBSplineCoefficientsToValuesGenericFilter_txx
  20
  21 /* =================================================
  22  * @file   clitkBSplineCoefficientsToValuesGenericFilter.txx
  23  * @author
  24  * @date
  25  *
  26  * @brief
  27  *
  28  ===================================================*/
  29
  30
  31 namespace clitk
  32 {
  33
  34   //-------------------------------------------------------------------
  35   // Update with the number of dimensions
  36   //-------------------------------------------------------------------
  37   template<unsigned int Dimension>
  38   void
  39   BSplineCoefficientsToValuesGenericFilter::UpdateWithDim(std::string PixelType, unsigned int Components)
  40   {
  41     if (m_Verbose) std::cout << "Image was detected to be "<<Dimension<<"D and "<<Components<<" component(s) of "<<  PixelType<<"..."<<std::endl;
  42
  43     if (Components==1)
  44       {
  45         //      if(PixelType == "short"){
  46         //        if (m_Verbose) std::cout << "Launching filter in "<< Dimension <<"D and signed short..." << std::endl;
  47         //        UpdateWithDimAndPixelType<Dimension, signed short>();
  48         //      }
  49         //    else if(PixelType == "unsigned_short"){
  50         //       if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and unsigned_short..." << std::endl;
  51         //       UpdateWithDimAndPixelType<Dimension, unsigned short>();
  52         //     }
  53
  54         //      else if (PixelType == "unsigned_char"){
  55         //        if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and unsigned_char..." << std::endl;
  56         //        UpdateWithDimAndPixelType<Dimension, unsigned char>();
  57         //      }
  58
  59         //     else if (PixelType == "char"){
  60         //       if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and signed_char..." << std::endl;
  61         //       UpdateWithDimAndPixelType<Dimension, signed char>();
  62         //     }
  63         //      else {
  64           if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and float..." << std::endl;
  65           UpdateWithDimAndPixelType<Dimension, float>();
  66           //    }
  67       }
  68
  69    else
  70     if (Components==3)
  71       {
  72         if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and 3D float (DVF)" << std::endl;
  73         UpdateWithDimAndVectorType<Dimension, itk::Vector<float, 3> >();
  74       }
  75
  76     else std::cerr<<"Number of components is "<<Components<<", not supported!"<<std::endl;
  77
  78   }
  79
  80
  81   //-------------------------------------------------------------------
  82   // Update with the number of dimensions and the pixeltype
  83   //-------------------------------------------------------------------
  84   template <unsigned int Dimension, class  PixelType>
  85   void
  86   BSplineCoefficientsToValuesGenericFilter::UpdateWithDimAndPixelType()
  87   {
  88
  89     // ImageTypes
  90     typedef itk::Image<PixelType, Dimension> InputImageType;
  91     typedef itk::Image<PixelType, Dimension> OutputImageType;
  92
  93     // Read the input
  94     typedef itk::ImageFileReader<InputImageType> InputReaderType;
  95     typename InputReaderType::Pointer reader = InputReaderType::New();
  96     reader->SetFileName( m_InputFileName);
  97     reader->Update();
  98     typename InputImageType::Pointer input= reader->GetOutput();
  99
 100     //Filter
 101     typedef  itk::ResampleImageFilter< InputImageType,OutputImageType >  ResampleFilterType;
 102     typename ResampleFilterType::Pointer resampler = ResampleFilterType::New();
 103
 104     // Properties
 105     if (m_ArgsInfo.like_given)
 106       {
 107         typename InputReaderType::Pointer likeReader=InputReaderType::New();
 108         likeReader->SetFileName(m_ArgsInfo.like_arg);
 109         likeReader->Update();
 110         resampler->SetOutputParametersFromImage(likeReader->GetOutput());
 111         if (m_Verbose)std::cout<<"Resampling output like "<<m_ArgsInfo.like_arg<<"..."<<std::endl;
 112       }
 113     else
 114       {
 115         // Size
 116         typename OutputImageType::SizeType outputSize;
 117         if (m_ArgsInfo.size_given)
 118           {
 119             for(unsigned int i=0; i< Dimension; i++)
 120               outputSize[i]=m_ArgsInfo.size_arg[i];
 121           }
 122         else outputSize=input->GetLargestPossibleRegion().GetSize();
 123         if (m_Verbose)std::cout<<"Setting the size to "<<outputSize<<"..."<<std::endl;
 124
 125         // Spacing
 126         typename OutputImageType::SpacingType outputSpacing;
 127         if (m_ArgsInfo.spacing_given)
 128           {
 129             for(unsigned int i=0; i< Dimension; i++)
 130               outputSpacing[i]=m_ArgsInfo.spacing_arg[i];
 131           }
 132         else outputSpacing=input->GetSpacing();
 133         if (m_Verbose)std::cout<<"Setting the spacing to "<<outputSpacing<<"..."<<std::endl;
 134
 135         // Origin
 136         typename OutputImageType::PointType outputOrigin;
 137         if (m_ArgsInfo.origin_given)
 138           {
 139             for(unsigned int i=0; i< Dimension; i++)
 140               outputOrigin[i]=m_ArgsInfo.origin_arg[i];
 141           }
 142         else outputOrigin=input->GetOrigin();
 143         if (m_Verbose)std::cout<<"Setting the origin to "<<outputOrigin<<"..."<<std::endl;
 144
 145         // Set
 146         resampler->SetSize( outputSize );
 147         resampler->SetOutputSpacing( outputSpacing );
 148         resampler->SetOutputOrigin(  outputOrigin );
 149
 150       }
 151
 152     // Interp : coeff
 153     typedef itk::BSplineResampleImageFunction<InputImageType, double> InterpolatorType;
 154     typename InterpolatorType::Pointer interpolator=InterpolatorType::New();
 155     interpolator->SetSplineOrder(m_ArgsInfo.order_arg);
 156
 157     // Set
 158     resampler->SetInterpolator( interpolator);
 159     resampler->SetDefaultPixelValue( static_cast<PixelType>(m_ArgsInfo.pad_arg) );
 160     resampler->SetInput( input );
 161
 162     // Update
 163     try
 164       {
 165         resampler->Update();
 166       }
 167     catch(itk::ExceptionObject)
 168       {
 169         std::cerr<<"Error resampling the image"<<std::endl;
 170       }
 171
 172     typename OutputImageType::Pointer output = resampler->GetOutput();
 173
 174     // Output
 175     typedef itk::ImageFileWriter<OutputImageType> WriterType;
 176     typename WriterType::Pointer writer = WriterType::New();
 177     writer->SetFileName(m_ArgsInfo.output_arg);
 178     writer->SetInput(output);
 179     writer->Update();
 180   }
 181
 182
 183   //-------------------------------------------------------------------
 184   // Update with the number of dimensions and the Vectortype
 185   //-------------------------------------------------------------------
 186   template <unsigned int Dimension, class  PixelType>
 187   void
 188   BSplineCoefficientsToValuesGenericFilter::UpdateWithDimAndVectorType()
 189   {
 190
 191     // ImageTypes
 192     typedef itk::Image<PixelType, Dimension> InputImageType;
 193     typedef itk::Image<PixelType, Dimension> OutputImageType;
 194
 195     // Read the input
 196     typedef itk::ImageFileReader<InputImageType> InputReaderType;
 197     typename InputReaderType::Pointer reader = InputReaderType::New();
 198     reader->SetFileName( m_InputFileName);
 199     reader->Update();
 200     typename InputImageType::Pointer input= reader->GetOutput();
 201
 202     //Filter
 203     typedef  itk::VectorResampleImageFilter< InputImageType,OutputImageType >  ResampleFilterType;
 204     typename ResampleFilterType::Pointer resampler = ResampleFilterType::New();
 205
 206     // Properties
 207     if (m_ArgsInfo.like_given)
 208       {
 209         typename InputReaderType::Pointer likeReader=InputReaderType::New();
 210         likeReader->SetFileName(m_ArgsInfo.like_arg);
 211         likeReader->Update();
 212
 213         // Set
 214         resampler->SetSize( likeReader->GetOutput()->GetLargestPossibleRegion().GetSize() );
 215         resampler->SetOutputSpacing(likeReader->GetOutput()->GetSpacing() );
 216         resampler->SetOutputOrigin(  likeReader->GetOutput()->GetOrigin() );
 217         if (m_Verbose)std::cout<<"Resampling output like "<<m_ArgsInfo.like_arg<<"..."<<std::endl;
 218       }
 219     else
 220       {
 221         // Size
 222         typename OutputImageType::SizeType outputSize;
 223         if (m_ArgsInfo.size_given)
 224           {
 225             for(unsigned int i=0; i< Dimension; i++)
 226               outputSize[i]=m_ArgsInfo.size_arg[i];
 227           }
 228         else outputSize=input->GetLargestPossibleRegion().GetSize();
 229         if (m_Verbose)std::cout<<"Setting the size to "<<outputSize<<"..."<<std::endl;
 230
 231         // Spacing
 232         typename OutputImageType::SpacingType outputSpacing;
 233         if (m_ArgsInfo.spacing_given)
 234           {
 235             for(unsigned int i=0; i< Dimension; i++)
 236               outputSpacing[i]=m_ArgsInfo.spacing_arg[i];
 237           }
 238         else outputSpacing=input->GetSpacing();
 239         if (m_Verbose)std::cout<<"Setting the spacing to "<<outputSpacing<<"..."<<std::endl;
 240
 241         // Origin
 242         typename OutputImageType::PointType outputOrigin;
 243         if (m_ArgsInfo.origin_given)
 244           {
 245             for(unsigned int i=0; i< Dimension; i++)
 246               outputOrigin[i]=m_ArgsInfo.origin_arg[i];
 247           }
 248         else outputOrigin=input->GetOrigin();
 249         if (m_Verbose)std::cout<<"Setting the origin to "<<outputOrigin<<"..."<<std::endl;
 250
 251         // Set
 252         resampler->SetSize( outputSize );
 253         resampler->SetOutputSpacing( outputSpacing );
 254         resampler->SetOutputOrigin(  outputOrigin );
 255
 256       }
 257
 258     // Interp : coeff
 259     typedef clitk::VectorBSplineResampleImageFunction<InputImageType, double> InterpolatorType;
 260     typename InterpolatorType::Pointer interpolator=InterpolatorType::New();
 261     interpolator->SetSplineOrder(m_ArgsInfo.order_arg);
 262
 263
 264     // Set
 265     resampler->SetInterpolator( interpolator);
 266     resampler->SetDefaultPixelValue( static_cast<PixelType>(m_ArgsInfo.pad_arg) );
 267     resampler->SetInput( input );
 268
 269     // Update
 270     try
 271       {
 272         resampler->Update();
 273       }
 274     catch(itk::ExceptionObject)
 275       {
 276         std::cerr<<"Error resampling the image"<<std::endl;
 277       }
 278
 279     typename OutputImageType::Pointer output = resampler->GetOutput();
 280     //writeImage<OutputImageType>(output, "/home/jef/tmp/output1.mhd", true);
 281
 282
 283     // Matrix Transform
 284     if(m_ArgsInfo.matrix_given)
 285       {
 286         typedef itk::TransformToDeformationFieldSource<OutputImageType, double> ConvertorType;
 287         typename   ConvertorType::Pointer filter= ConvertorType::New();
 288         filter->SetOutputParametersFromImage(output);
 289
 290         typedef itk::AffineTransform<double, Dimension> TransformType;
 291         typename TransformType::Pointer transform =TransformType::New();
 292         itk::Matrix<double, Dimension+1, Dimension+1> homMatrix= ReadMatrix<Dimension>( m_ArgsInfo.matrix_arg);
 293         itk::Matrix<double, Dimension, Dimension> matrix =GetRotationalPartMatrix( homMatrix);
 294         itk::Vector<double, Dimension> offset= GetTranslationPartMatrix( homMatrix);
 295         transform->SetMatrix(matrix);
 296         transform->SetOffset(offset);
 297         filter->SetTransform(transform);
 298         filter->Update();
 299         typename OutputImageType::Pointer output2=filter->GetOutput();
 300         //writeImage<OutputImageType>(output2, "/home/jef/tmp/output2.mhd", true);
 301
 302         // Add
 303         typedef itk::AddImageFilter< OutputImageType, OutputImageType, OutputImageType > AddType;
 304         typename AddType::Pointer adder= AddType::New();
 305         adder->SetInput(output);
 306         adder->SetInput(1,output2);
 307         adder->Update();
 308         output=adder->GetOutput();
 309       }
 310
 311
 312     // Output
 313     typedef itk::ImageFileWriter<OutputImageType> WriterType;
 314     typename WriterType::Pointer writer = WriterType::New();
 315     writer->SetFileName(m_ArgsInfo.output_arg);
 316     writer->SetInput(output);
 317     writer->Update();
 318   }
 319
 320 }//end clitk
 321
 322 #endif //#define clitkBSplineCoefficientsToValuesGenericFilter_txx