tools/clitkWriteDicomSeriesGenericFilter.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://www.centreleonberard.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 clitkWriteDicomSeriesGenericFilter_txx
  19 #define clitkWriteDicomSeriesGenericFilter_txx
  20
  21 /* =================================================
  22  * @file   clitkWriteDicomSeriesGenericFilter.txx
  23  * @author
  24  * @date
  25  *
  26  * @brief
  27  *
  28  ===================================================*/
  29
  30 // clitk
  31 #include "clitkResampleImageWithOptionsFilter.h"
  32
  33
  34 namespace clitk
  35 {
  36
  37
  38 //-----------------------------------------------------------
  39 // Constructor
  40 //-----------------------------------------------------------
  41 template<class args_info_type>
  42 WriteDicomSeriesGenericFilter<args_info_type>::WriteDicomSeriesGenericFilter()
  43 {
  44   m_Verbose=false;
  45   m_InputFileName="";
  46 }
  47
  48
  49 //-----------------------------------------------------------
  50 // Update
  51 //-----------------------------------------------------------
  52 template<class args_info_type>
  53 void WriteDicomSeriesGenericFilter<args_info_type>::Update()
  54 {
  55   // Read the Dimension and PixelType
  56   int Dimension;
  57   std::string PixelType;
  58   ReadImageDimensionAndPixelType(m_InputFileName, Dimension, PixelType);
  59
  60
  61   // Call UpdateWithDim
  62   if(Dimension==2) UpdateWithDim<2>(PixelType);
  63   else if(Dimension==3) UpdateWithDim<3>(PixelType);
  64   // else if (Dimension==4)UpdateWithDim<4>(PixelType);
  65   else {
  66     std::cout<<"Error, Only for 2 or 3  Dimensions!!!"<<std::endl ;
  67     return;
  68   }
  69 }
  70
  71 //-------------------------------------------------------------------
  72 // Update with the number of dimensions
  73 //-------------------------------------------------------------------
  74 template<class args_info_type>
  75 template<unsigned int Dimension>
  76 void
  77 WriteDicomSeriesGenericFilter<args_info_type>::UpdateWithDim(std::string PixelType)
  78 {
  79   if (m_Verbose) std::cout << "Image was detected to be "<<Dimension<<"D and "<< PixelType<<"..."<<std::endl;
  80
  81   if(PixelType == "short") {
  82     if (m_Verbose) std::cout << "Launching filter in "<< Dimension <<"D and signed short..." << std::endl;
  83     UpdateWithDimAndPixelType<Dimension, signed short>();
  84   }
  85   //    else if(PixelType == "unsigned_short"){
  86   //       if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and unsigned_short..." << std::endl;
  87   //       UpdateWithDimAndPixelType<Dimension, unsigned short>();
  88   //     }
  89
  90   else if (PixelType == "unsigned_char") {
  91     if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and unsigned_char..." << std::endl;
  92     UpdateWithDimAndPixelType<Dimension, unsigned char>();
  93   }
  94
  95   //     else if (PixelType == "char"){
  96   //       if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and signed_char..." << std::endl;
  97   //       UpdateWithDimAndPixelType<Dimension, signed char>();
  98   //     }
  99   else {
 100     if (m_Verbose) std::cout  << "Launching filter in "<< Dimension <<"D and float..." << std::endl;
 101     UpdateWithDimAndPixelType<Dimension, float>();
 102   }
 103 }
 104
 105 //-------------------------------------------------------------------
 106 // Update with the number of dimensions and the pixeltype read from
 107 // the dicom files. The MHD files may be resampled to match the
 108 // dicom spacing (and number of slices). Rounding errors in resampling
 109 // are handled by removing files when generating the output dicom
 110 // series.
 111 //-------------------------------------------------------------------
 112 template<class args_info_type>
 113 template <unsigned int Dimension, class  PixelType>
 114 void
 115 WriteDicomSeriesGenericFilter<args_info_type>::UpdateWithDimAndPixelType()
 116 {
 117
 118   // ImageTypes
 119   typedef itk::Image<PixelType, Dimension> InputImageType;
 120   typedef itk::Image<PixelType, 2> OutputImageType;
 121
 122   // Read the dicom directory
 123   typedef itk::ImageSeriesReader< InputImageType >     ReaderType;
 124   typedef itk::GDCMImageIO ImageIOType;
 125   typedef itk::GDCMSeriesFileNames NamesGeneratorType;
 126
 127   ImageIOType::Pointer gdcmIO = ImageIOType::New();
 128   NamesGeneratorType::Pointer namesGenerator = NamesGeneratorType::New();
 129   namesGenerator->SetInputDirectory( m_ArgsInfo.inputDir_arg );
 130   namesGenerator->SetOutputDirectory( m_ArgsInfo.outputDir_arg  );
 131   typename   ReaderType::FileNamesContainer filenames_in = namesGenerator->GetInputFileNames();
 132   typename   ReaderType::FileNamesContainer filenames_out = namesGenerator->GetOutputFileNames();
 133
 134   // Output the dicom files
 135   unsigned int numberOfFilenames =  filenames_in.size();
 136   if (m_Verbose) {
 137     std::cout << numberOfFilenames <<" were read in the directory "<<m_ArgsInfo.inputDir_arg<<"..."<<std::endl<<std::endl;
 138     for(unsigned int fni = 0; fni<numberOfFilenames; fni++) {
 139       std::cout << "filename # " << fni << " = ";
 140       std::cout << filenames_in[fni] << std::endl;
 141     }
 142   }
 143
 144   // Read the series
 145   typename ReaderType::Pointer reader = ReaderType::New();
 146   reader->SetImageIO( gdcmIO );
 147   reader->SetFileNames( filenames_in );
 148   try {
 149     reader->Update();
 150   } catch (itk::ExceptionObject &excp) {
 151     std::cerr << "Error: Exception thrown while reading the DICOM series!!" << std::endl;
 152     std::cerr << excp << std::endl;
 153   }
 154
 155   // Read the input (MHD file)
 156   typedef typename InputImageType::RegionType RegionType;
 157   typedef typename RegionType::SizeType SizeType;
 158   typedef itk::ImageFileReader<InputImageType> InputReaderType;
 159   typename InputReaderType::Pointer volumeReader = InputReaderType::New();
 160   volumeReader->SetFileName( m_InputFileName);
 161   volumeReader->Update();
 162
 163   typename InputImageType::Pointer input = volumeReader->GetOutput();
 164   if ((!m_ArgsInfo.useSizeAsReference_flag && (input->GetSpacing() != reader->GetOutput()->GetSpacing())) ||
 165       (m_ArgsInfo.useSizeAsReference_flag && (input->GetLargestPossibleRegion().GetSize() != reader->GetOutput()->GetLargestPossibleRegion().GetSize()))) {
 166
 167     // resampling is carried out on the fly if resolution or size between
 168     // the input mhd and input dicom series is different
 169
 170     // Filter
 171     typedef clitk::ResampleImageWithOptionsFilter<InputImageType, InputImageType> ResampleImageFilterType;
 172     typename ResampleImageFilterType::Pointer filter = ResampleImageFilterType::New();
 173     filter->SetInput(input);
 174     filter->SetVerboseOptions(m_Verbose);
 175     filter->SetGaussianFilteringEnabled(false);
 176     filter->SetDefaultPixelValue(0);
 177
 178     if (!m_ArgsInfo.useSizeAsReference_flag) {
 179       filter->SetOutputSpacing(reader->GetOutput()->GetSpacing());
 180       if (m_Verbose) {
 181         std::cout << "Warning: The image spacing differs between the MHD file and the input dicom series. Performing resampling with default options using spacing as reference (for advanced options, use clitkResampleImage)." << std::endl;
 182         std::cout << "MHD -> " << input->GetSpacing() << std::endl;
 183         std::cout << "dicom -> " << reader->GetOutput()->GetSpacing() << std::endl;
 184       }
 185     }
 186     else {
 187       const SizeType& dicom_size = reader->GetOutput()->GetLargestPossibleRegion().GetSize();
 188       SizeType output_size;
 189       for (unsigned int i = 0; i < Dimension; i++)
 190         output_size[i] = dicom_size[i];
 191       filter->SetOutputSize(output_size);
 192       if (m_Verbose) {
 193           std::cout << "Warning: The image size differs between the MHD file and the input dicom series. Performing resampling with default options using size as reference (for advanced options, use clitkResampleImage)." << std::endl;
 194           std::cout << "MHD -> " << input->GetLargestPossibleRegion().GetSize() << std::endl;
 195           std::cout << "dicom -> " << reader->GetOutput()->GetLargestPossibleRegion().GetSize() << std::endl;
 196       }
 197     }
 198
 199     try {
 200       filter->Update();
 201       input = filter->GetOutput();
 202     } catch( itk::ExceptionObject & excp ) {
 203     std::cerr << "Error: Exception thrown while resampling!!" << std::endl;
 204     std::cerr << excp << std::endl;
 205     }
 206   }
 207
 208   //    In some cases, due to resampling approximation issues,
 209   //    the number of slices in the MHD file may be different (smaller)
 210   //    from the number of files in the template dicom directory.
 211   //    To avoid ITK generating an exception, we reduce the number
 212   //    of DCM files to be considered, and a warning is printed
 213   //    in verbose mode
 214   const RegionType volumeRegion = input->GetLargestPossibleRegion();
 215   const SizeType& volumeSize = volumeRegion.GetSize();
 216   if (m_Verbose) {
 217     std::cout << volumeRegion << volumeSize << std::endl;
 218   }
 219   if (Dimension == 3 && volumeSize[2] < numberOfFilenames) {
 220     if (m_Verbose)
 221       std::cout << "Warning: The number of files in " << m_ArgsInfo.inputDir_arg << " (" << filenames_in.size() << " files) is greater than the number of slices in MHD (" << volumeSize[2] << " slices). Using only " << volumeSize[2] << " files." << std::endl;
 222
 223     filenames_in.resize(volumeSize[2]);
 224     filenames_out.resize(filenames_in.size());
 225     numberOfFilenames =  filenames_in.size();
 226   }
 227
 228   // Modify the meta dictionary
 229   typedef itk::MetaDataDictionary   DictionaryType;
 230   const std::vector<DictionaryType*>* dictionary = reader->GetMetaDataDictionaryArray();
 231
 232   // Get keys
 233   unsigned int numberOfKeysGiven=0;
 234   if(m_ArgsInfo.midP_flag && m_ArgsInfo.key_given)
 235     std::cerr<<"Error: both keys and midP option are given"<<std::endl;
 236   else if (m_ArgsInfo.midP_flag)
 237     numberOfKeysGiven=1;
 238   else
 239     numberOfKeysGiven=m_ArgsInfo.key_given;
 240
 241   for (unsigned int i = 0; i < numberOfKeysGiven; i++) {
 242     std::string entryId(m_ArgsInfo.key_arg[i]  );
 243     std::string value( m_ArgsInfo.tag_arg[i] );
 244     std::cout << entryId << " " << value << std::endl;
 245     for(unsigned int fni = 0; fni<numberOfFilenames; fni++)
 246       itk::EncapsulateMetaData<std::string>( *((*dictionary)[fni]), entryId, value );
 247   }
 248
 249
 250   // Output directory and filenames
 251   itksys::SystemTools::MakeDirectory( m_ArgsInfo.outputDir_arg ); // create if it doesn't exist
 252   typedef itk::ImageSeriesWriter<InputImageType, OutputImageType >  SeriesWriterType;
 253   typename SeriesWriterType::Pointer seriesWriter = SeriesWriterType::New();
 254
 255   seriesWriter->SetInput( input );
 256   seriesWriter->SetImageIO( gdcmIO );
 257
 258   seriesWriter->SetFileNames( filenames_out );
 259   seriesWriter->SetMetaDataDictionaryArray( dictionary );
 260
 261   // Write
 262   try {
 263     seriesWriter->Update();
 264     std::cerr << "filter update" << std::endl;
 265   } catch( itk::ExceptionObject & excp ) {
 266     std::cerr << "Error: Exception thrown while writing the series!!" << std::endl;
 267     std::cerr << excp << std::endl;
 268   }
 269
 270 }
 271
 272 /*
 273 //-------------------------------------------------------------------
 274 // Update with the number of dimensions and the pixeltype
 275 //-------------------------------------------------------------------
 276 template<class args_info_type>
 277 template <unsigned int Dimension, class  PixelType>
 278 void
 279 WriteDicomSeriesGenericFilter<args_info_type>::UpdateWithDimAndPixelType()
 280 {
 281
 282   // ImageTypes
 283   typedef itk::Image<PixelType, Dimension> InputImageType;
 284   typedef itk::Image<PixelType, 2> OutputImageType;
 285
 286   // Read the input (volumetric)
 287   typedef itk::ImageFileReader<InputImageType> InputReaderType;
 288   typename InputReaderType::Pointer volumeReader = InputReaderType::New();
 289   volumeReader->SetFileName( m_InputFileName);
 290   volumeReader->Update();
 291   typename InputImageType::Pointer input= volumeReader->GetOutput();
 292
 293   // Read the dicom directory
 294   typedef itk::ImageSeriesReader< InputImageType >     ReaderType;
 295   typedef itk::GDCMImageIO ImageIOType;
 296   typedef itk::GDCMSeriesFileNames NamesGeneratorType;
 297
 298   ImageIOType::Pointer gdcmIO = ImageIOType::New();
 299   NamesGeneratorType::Pointer namesGenerator = NamesGeneratorType::New();
 300   namesGenerator->SetInputDirectory( m_ArgsInfo.inputDir_arg );
 301   namesGenerator->SetOutputDirectory( m_ArgsInfo.outputDir_arg  );
 302   typename   ReaderType::FileNamesContainer filenames_in = namesGenerator->GetInputFileNames();
 303   typename   ReaderType::FileNamesContainer filenames_out = namesGenerator->GetOutputFileNames();
 304
 305   // Output the dicom files
 306   unsigned int numberOfFilenames =  filenames_in.size();
 307   if (m_Verbose) {
 308     std::cout << numberOfFilenames <<" were read in the directory "<<m_ArgsInfo.inputDir_arg<<"..."<<std::endl<<std::endl;
 309     for(unsigned int fni = 0; fni<numberOfFilenames; fni++) {
 310       std::cout << "filename # " << fni << " = ";
 311       std::cout << filenames_in[fni] << std::endl;
 312     }
 313   }
 314
 315   // RP: 16/03/2011
 316   //    In some cases, due to resampling approximation issues,
 317   //    the number of slices in the MHD file may be different
 318   //    from the number of slices in the template DCM directory.
 319   //    To avoid ITK generating an exception, we reduce the number
 320   //    of DCM files to be considered, provided the --force
 321   //    option is set.
 322   typedef typename InputImageType::RegionType RegionType;
 323   typedef typename RegionType::SizeType SizeType;
 324   const RegionType volumeRegion = input->GetLargestPossibleRegion();
 325   const SizeType& volumeSize = volumeRegion.GetSize();
 326   if (m_ArgsInfo.force_given && Dimension == 3 && volumeSize[2] < numberOfFilenames)
 327   {
 328     std::cout << "Warning: Number of files in " << m_ArgsInfo.inputDir_arg << " is greater than the number of slices in MHD. Using only " << volumeSize[2] << " files." << std::endl;
 329     filenames_in.resize(volumeSize[2]);
 330     filenames_out.resize(filenames_in.size());
 331     numberOfFilenames =  filenames_in.size();
 332   }
 333
 334   // Read the series
 335   typename ReaderType::Pointer reader = ReaderType::New();
 336   reader->SetImageIO( gdcmIO );
 337   reader->SetFileNames( filenames_in );
 338   try {
 339     reader->Update();
 340   } catch (itk::ExceptionObject &excp) {
 341     std::cerr << "Error: Exception thrown while reading the DICOM series!!" << std::endl;
 342     std::cerr << excp << std::endl;
 343   }
 344
 345
 346   // Modify the meta dictionary
 347   typedef itk::MetaDataDictionary   DictionaryType;
 348   const std::vector<DictionaryType*>* dictionary = reader->GetMetaDataDictionaryArray();
 349
 350   // Get keys
 351   unsigned int numberOfKeysGiven=0;
 352   if(m_ArgsInfo.midP_flag && m_ArgsInfo.key_given)
 353     std::cerr<<"Error: both keys and midP option are given"<<std::endl;
 354   else if (m_ArgsInfo.midP_flag)
 355     numberOfKeysGiven=1;
 356   else
 357     numberOfKeysGiven=m_ArgsInfo.key_given;
 358
 359   for (unsigned int i = 0; i < numberOfKeysGiven; i++) {
 360     std::string entryId(m_ArgsInfo.key_arg[i]  );
 361     std::string value( m_ArgsInfo.tag_arg[i] );
 362     for(unsigned int fni = 0; fni<numberOfFilenames; fni++)
 363       itk::EncapsulateMetaData<std::string>( *((*dictionary)[fni]), entryId, value );
 364   }
 365
 366   // Output directory and filenames
 367   itksys::SystemTools::MakeDirectory( m_ArgsInfo.outputDir_arg ); // create if it doesn't exist
 368   typedef itk::ImageSeriesWriter<InputImageType, OutputImageType >  SeriesWriterType;
 369   typename SeriesWriterType::Pointer seriesWriter = SeriesWriterType::New();
 370
 371   seriesWriter->SetInput( volumeReader->GetOutput() );
 372   seriesWriter->SetImageIO( gdcmIO );
 373
 374   seriesWriter->SetFileNames( filenames_out );
 375   seriesWriter->SetMetaDataDictionaryArray( dictionary );
 376
 377   // Write
 378   try {
 379     seriesWriter->Update();
 380   } catch( itk::ExceptionObject & excp ) {
 381     std::cerr << "Error: Exception thrown while writing the series!!" << std::endl;
 382     std::cerr << excp << std::endl;
 383   }
 384
 385 }
 386 */
 387
 388 }//end clitk
 389
 390 #endif //#define clitkWriteDicomSeriesGenericFilter_txx