1 /*=========================================================================
2 Program: vv http://www.creatis.insa-lyon.fr/rio/vv
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
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.
13 It is distributed under dual licence
15 - BSD See included LICENSE.txt file
16 - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html
17 ===========================================================================**/
18 #ifndef clitkImage2DicomDoseGenericFilter_txx
19 #define clitkImage2DicomDoseGenericFilter_txx
21 /* =================================================
22 * @file clitkImage2DicomDosemGenericFilter.txx
28 ===================================================*/
32 #include "clitkImage2DicomDoseGenericFilter.h"
33 #include "clitkCommon.h"
34 #if GDCM_MAJOR_VERSION >= 2
35 #include "gdcmUIDGenerator.h"
36 #include <gdcmImageHelper.h>
37 #include <gdcmAttribute.h>
38 #include <gdcmReader.h>
39 #include <gdcmWriter.h>
45 //#include "gdcmBase.h"
46 //#include "gdcmDocEntry.h"
54 //-----------------------------------------------------------
56 //-----------------------------------------------------------
57 template<class args_info_type>
58 Image2DicomDoseGenericFilter<args_info_type>::Image2DicomDoseGenericFilter()
65 //-----------------------------------------------------------
67 //-----------------------------------------------------------
68 template<class args_info_type>
69 void Image2DicomDoseGenericFilter<args_info_type>::Update()
71 // Read the Dimension and PixelType
73 std::string PixelType;
74 ReadImageDimensionAndPixelType(m_InputFileName, Dimension, PixelType);
78 if(Dimension==2) UpdateWithDim<2>(PixelType);
79 else if(Dimension==3) UpdateWithDim<3>(PixelType);
80 // else if (Dimension==4)UpdateWithDim<4>(PixelType);
82 std::cout<<"Error, Only for 2 or 3 Dimensions!!!"<<std::endl ;
87 //-------------------------------------------------------------------
88 // Update with the number of dimensions
89 //-------------------------------------------------------------------
90 template<class args_info_type>
91 template<unsigned int Dimension>
93 Image2DicomDoseGenericFilter<args_info_type>::UpdateWithDim(std::string PixelType)
95 if (m_Verbose) std::cout << "Image was detected to be "<<Dimension<<"D and "<< PixelType<<"..."<<std::endl;
97 if(PixelType == "short") {
98 if (m_Verbose) std::cout << "Launching filter in "<< Dimension <<"D and signed short..." << std::endl;
99 UpdateWithDimAndPixelType<Dimension, signed short>();
101 else if(PixelType == "unsigned_short"){
102 if (m_Verbose) std::cout << "Launching filter in "<< Dimension <<"D and unsigned_short..." << std::endl;
103 UpdateWithDimAndPixelType<Dimension, unsigned short>();
106 else if (PixelType == "unsigned_char") {
107 if (m_Verbose) std::cout << "Launching filter in "<< Dimension <<"D and unsigned_char..." << std::endl;
108 UpdateWithDimAndPixelType<Dimension, unsigned char>();
111 // else if (PixelType == "char"){
112 // if (m_Verbose) std::cout << "Launching filter in "<< Dimension <<"D and signed_char..." << std::endl;
113 // UpdateWithDimAndPixelType<Dimension, signed char>();
115 else if (PixelType == "double") {
116 if (m_Verbose) std::cout << "Launching filter in "<< Dimension <<"D and double..." << std::endl;
117 UpdateWithDimAndPixelType<Dimension, double>();
120 if (m_Verbose) std::cout << "Launching filter in "<< Dimension <<"D and float..." << std::endl;
121 UpdateWithDimAndPixelType<Dimension, float>();
125 //-------------------------------------------------------------------
126 // Update with the number of dimensions and the pixeltype read from
127 // the dicom files. The MHD files may be resampled to match the
128 // dicom spacing (and number of slices). Rounding errors in resampling
129 // are handled by removing files when generating the output dicom
131 //-------------------------------------------------------------------
132 template<class args_info_type>
133 template <unsigned int Dimension, class PixelType>
135 Image2DicomDoseGenericFilter<args_info_type>::UpdateWithDimAndPixelType()
138 #if GDCM_MAJOR_VERSION == 2
140 typedef itk::Image<PixelType, Dimension> InputImageType;
141 typedef itk::Image<PixelType, Dimension> OutputImageType;
142 typedef itk::ImageFileReader< InputImageType > ReaderType;
143 typedef itk::ImageSeriesReader< InputImageType > ReaderSeriesType;
144 typedef itk::ImageRegionIterator< InputImageType > IteratorType;
145 typedef itk::GDCMImageIO ImageIOType;
147 //-----------------------------------------------------------------------------
148 // opening image input file
149 typename ReaderType::Pointer reader = ReaderType::New();
150 const char * filename = m_ArgsInfo.input_arg;
151 reader->SetFileName( filename );
153 typename InputImageType::Pointer image = reader->GetOutput();
156 typename InputImageType::PointType origin = image->GetOrigin();
160 typename InputImageType::SizeType imageSize = image->GetLargestPossibleRegion().GetSize();
162 int NbCols=imageSize[0]; // col
163 int NbRows=imageSize[1]; // row
164 int NbFrames=imageSize[2]; // frame
170 typename InputImageType::SpacingType Spacing = image->GetSpacing();
175 std::stringstream strOffset;
177 for (int i=1; i<NbFrames ; i++){
182 DD(strOffset.str().c_str());
185 float highestValue=pow(10,-10);
186 IteratorType out( image, image->GetRequestedRegion() );
187 for (out.GoToBegin(); !out.IsAtEnd(); ++out){
189 if (out.Get()>highestValue) highestValue=out.Get();
191 double doseScaling = highestValue/(pow(2,16)-1);
195 /*std::vector<unsigned short int> ImageData;
196 typename InputImageType::IndexType pixelIndex;
198 unsigned short int pixelValue;
200 for (int i=0; i<NbFrames; i++){
202 for (int j=0; j<NbRows; j++){
204 for (int k=0; k<NbCols; k++){
206 pixelValue=image->GetPixel(pixelIndex)/doseScaling;
207 if(float(image->GetPixel(pixelIndex)/doseScaling)>(pow(2,16)-1.)) {
208 std::cout<<"\n!!!!! WARNING !!!!! pixel index: "<<pixelIndex<<"unsigned short int capacity ful or overfuled => Highest value may become 0"<<std::endl;
210 DD(image->GetPixel(pixelIndex));
211 //DD(image->GetPixel(pixelIndex)/doseScaling);
213 std::cout<<"Pixel Value should be equal to "<<(pow(2,16)-1)<<" but should not be 0"<<std::endl;
214 std::cout<<"\n"<<std::endl;
215 //assert(pixelValue<=(pow(2,16)-1)); should work, but do not...
218 ImageData.push_back(pixelValue);
223 DD(ImageData.size());*/
225 // Read Dicom model file
226 typename ReaderSeriesType::Pointer readerSeries = ReaderSeriesType::New();
227 ImageIOType::Pointer gdcmIO = ImageIOType::New();
228 std::string filename_out = m_ArgsInfo.output_arg;
229 gdcmIO->LoadPrivateTagsOn();
230 gdcmIO->KeepOriginalUIDOn();
231 typename ReaderSeriesType::FileNamesContainer fileNames;
232 fileNames.push_back(m_ArgsInfo.DicomInputFile_arg);
233 readerSeries->SetImageIO( gdcmIO );
234 readerSeries->SetFileNames( fileNames );
236 readerSeries->Update();
237 } catch (itk::ExceptionObject &excp) {
238 std::cerr << "Error: Exception thrown while reading the DICOM model file !!" << std::endl;
239 std::cerr << excp << std::endl;
242 // update output dicom keys/tags
243 // string for distinguishing items inside sequence:
244 const std::string ITEM_ENCAPSULATE_STRING("DICOM_ITEM_ENCAPSULATE");
245 std::string tempString = ITEM_ENCAPSULATE_STRING + "01";
246 typename ReaderSeriesType::DictionaryRawPointer inputDict = (*(readerSeries->GetMetaDataDictionaryArray()))[0];
247 typename ReaderSeriesType::DictionaryArrayType outputArray;
248 typename ReaderSeriesType::DictionaryRawPointer outputDict = new typename ReaderSeriesType::DictionaryType;
249 CopyDictionary (*inputDict, *outputDict);
250 itk::EncapsulateMetaData<std::string>(*outputDict, "0020|0032", NumberToString(origin));
251 itk::EncapsulateMetaData<std::string>(*outputDict, "0028|0008", NumberToString(NbFrames));
252 itk::EncapsulateMetaData<std::string>(*outputDict, "0028|0010", NumberToString(NbRows));
253 itk::EncapsulateMetaData<std::string>(*outputDict, "0028|0011", NumberToString(NbCols));
254 itk::EncapsulateMetaData<std::string>(*outputDict, "0028|0030", NumberToString(Spacing));
255 itk::EncapsulateMetaData<std::string>(*outputDict, "3004|000e", NumberToString(doseScaling));
256 itk::EncapsulateMetaData<std::string>(*outputDict, "3004|000c", strOffset.str());
257 outputArray.push_back(outputDict);
259 // Output directory and filenames
260 typedef itk::ImageSeriesWriter<OutputImageType, OutputImageType> WriterSerieType;
261 typename WriterSerieType::Pointer writerSerie = WriterSerieType::New();
262 writerSerie->SetInput( image );
263 writerSerie->SetImageIO( gdcmIO );
264 typename ReaderSeriesType::FileNamesContainer fileNamesOutput;
265 fileNamesOutput.push_back(filename_out);
266 writerSerie->SetFileNames( fileNamesOutput );
267 writerSerie->SetMetaDataDictionaryArray(&outputArray);
271 if (m_ArgsInfo.verbose_flag)
272 std::cout << writerSerie << std::endl;
273 writerSerie->Update();
274 } catch( itk::ExceptionObject & excp ) {
275 std::cerr << "Error: Exception thrown while writing the series!!" << std::endl;
276 std::cerr << excp << std::endl;
280 //---------------------------------------------------------------------------------------
282 // The previous way of writting DICOM-RT-DOSE works only for ITK
283 // and resulting RT-DOSE files are not readable by commercial systems.
284 // The nex step is to copy again the output file with another syntax, allowing to make RT-DOSE files readable by commercial systems.
285 // see Jean-Pierre Roux comment below.
287 /*gdcm::FileHelper *fh = new gdcm::FileHelper(args_info.OutputFile_arg);
291 dataSize = fh->GetImageDataRawSize();
292 imageData = fh->GetImageDataRaw();// somewhat important : Loads the Pixels in memory !
294 fh->SetWriteModeToRaw();
295 fh->SetWriteTypeToDcmExplVR();
297 bool res = fh->Write(args_info.OutputFile_arg);
300 std::cout <<"Fail to write [" << args_info.OutputFile_arg << "]" <<std::endl;
301 else std::cout<<"\n DICOM File re-written, using the FileHelper syntax, in order to be processed by commercial systems !"<<std::endl;
306 w.SetFileName(m_ArgsInfo.output_arg);
310 //---------------------------------------------------------------------------------------
311 //---------------------------------------------------------------------------------------
312 // Jean-Pierre Roux help for DICOM-RT-DOSE writting
313 //---------------------------------------------------------------------------------------
314 //---------------------------------------------------------------------------------------
316 // de Jean-Pierre Roux <jpr@creatis.insa-lyon.fr>
317 // répondre à jpr@creatis.insa-lyon.fr
318 // à Loic Grevillot <loic.grevillot@gmail.com>
319 // cc jpr@creatis.insa-lyon.fr,
320 // Joël Schaerer <joel.schaerer@gmail.com>,
321 // David Sarrut <David.Sarrut@creatis.insa-lyon.fr>,
322 // Joel Schaerer <joel.schaerer@creatis.insa-lyon.fr>
323 // date 12 juillet 2010 12:23
324 // objet Re: DICOM RT DOSE
326 // masquer les détails 12:23 (Il y a 21 heures)
330 // J'aurais écrit à peut prèt la même chose.
331 // (Ci après un extrait de mon code -Example/ReWrite.cxx-)
333 // L'utilisation d'un FileHelper (qui ne changera rien dans ce cas précis) est une mesure de précaution, car, l'élément 7FE0|0010 peut être compressé (ce qui n'est pas le cas pour tes images), que la manière de stocker les pixels ainsi compressés était parfois un peu ... curieuse.
334 // Je décompresse, et réécrit non compressé.
336 // ============================
337 // gdcm::File *f = new gdcm::File();
338 // f->SetMaxSizeLoadEntry(0x7fffffff);
339 // f->SetFileName( fileName );
340 // bool res = f->Load();
347 // if (!f->IsReadable())
349 // std::cerr << "Sorry, not a Readable DICOM / ACR File" <<std::endl;
354 // gdcm::FileHelper *fh = new gdcm::FileHelper(f);
358 // dataSize = fh->GetImageDataRawSize();
359 // imageData = fh->GetImageDataRaw();// somewhat important : Loads the Pixels in memory !
361 // fh->SetWriteModeToRaw();
362 // fh->SetWriteTypeToDcmExplVR();
363 // res = fh->Write(outputFileName);
366 // std::cout <<"Fail to write [" << outputFileName << "]" <<std::endl;
371 //---------------------------------------------------------------------------------------
372 //---------------------------------------------------------------------------------------
373 // pour supprimer des tags:
374 //---------------------------------------------------------------------------------------
375 //---------------------------------------------------------------------------------------
376 // SOLUCE 1 QUI MARCHE POUR UN SQItem SIMPLE
379 //a= ((gdcm::SQItem*)mDCMFile)->GetDocEntry(0xfffe,0xe00d);
380 a= ((gdcm::SQItem*)mDCMFile)->GetDocEntry(0x3004,0x000e);
381 ((gdcm::SQItem*)mDCMFile)->RemoveEntry(a);
382 mDCMFile->Write (args_info.OutputFile_arg, type);
385 // SOLUCE 2 QUI MARCHE POUR UNE SeqEntry->SQItem
387 std::cout<<"\ntest correction fichier apres ecriture\n"<<std::endl;
388 gdcm::SeqEntry *seqEntry = mDCMFile->GetSeqEntry(0x300c,0x0002);
389 gdcm::SQItem* currentItem = seqEntry->GetFirstSQItem();
391 //a= currentItem->GetDocEntry(0x0008,0x1155);
392 a= currentItem->GetDocEntry(0xfffe,0xe00d);
393 currentItem->RemoveEntry(a);
394 mDCMFile->Write (args_info.OutputFile_arg, type);
398 //a=GetDocEntry(0x7fe0,0x0000);
399 //((gdcm::SQItem*)mDCMFile)->RemoveEntry(a);
401 //-----------------------------------------------------------------------------------
402 std::cout <<"\n## DICOM Image to RT DOSE application is ended..."<<std::endl;
403 std::cout <<"#########################################################\n" << std::endl;
406 std::cout << "Use GDCM2" << std::endl;
409 //---------------------------------------------------------------------------
412 //---------------------------------------------------------------------------
413 void CopyDictionary (itk::MetaDataDictionary &fromDict, itk::MetaDataDictionary &toDict)
415 typedef itk::MetaDataDictionary DictionaryType;
417 DictionaryType::ConstIterator itr = fromDict.Begin();
418 DictionaryType::ConstIterator end = fromDict.End();
419 typedef itk::MetaDataObject< std::string > MetaDataStringType;
423 itk::MetaDataObjectBase::Pointer entry = itr->second;
425 MetaDataStringType::Pointer entryvalue =
426 dynamic_cast<MetaDataStringType *>( entry.GetPointer() ) ;
429 std::string tagkey = itr->first;
430 std::string tagvalue = entryvalue->GetMetaDataObjectValue();
431 itk::EncapsulateMetaData<std::string>(toDict, tagkey, tagvalue);
436 //---------------------------------------------------------------------------
438 //---------------------------------------------------------------------------
439 template <typename T> std::string NumberToString ( T Number )
441 std::ostringstream ss;
445 //---------------------------------------------------------------------------
450 #endif //#define clitkImage2DicomDoseGenericFilter_txx