#include <gdcmAttribute.h>
#include <gdcmReader.h>
#include <gdcmWriter.h>
+#include <gdcmDataElement.h>
+#include <gdcmTag.h>
#else
#include "gdcmFile.h"
#include "gdcmUtil.h"
// Read Dicom model file
typename ReaderType::Pointer reader = ReaderType::New();
ImageIOType::Pointer gdcmIO = ImageIOType::New();
- std::string filename_out = m_ArgsInfo.output_arg;
gdcmIO->LoadPrivateTagsOn();
gdcmIO->KeepOriginalUIDOn();
reader->SetImageIO( gdcmIO );
typename ReaderType::DictionaryArrayType outputArray;
typename ReaderType::DictionaryRawPointer outputDict = new typename ReaderType::DictionaryType;
CopyDictionary (*inputDict, *outputDict);
- //itk::EncapsulateMetaData<std::string>(*outputDict, "0054|0011", NumberToString(energyNumber));
- //itk::EncapsulateMetaData<std::string>(*outputDict, "0054|0021", NumberToString(headNumber));
+ itk::EncapsulateMetaData<std::string>(*outputDict, "0054|0011", NumberToString(1));
outputArray.push_back(outputDict);
- //Read the number of slices and energy window
+ //Read the number of slices and energy window, rotation
int nbSlice = input->GetLargestPossibleRegion().GetSize()[2];
gdcm::Reader hreader;
hreader.SetFileName(m_ArgsInfo.input_arg);
hreader.Read();
- gdcm::DataSet& ds = hreader.GetFile().GetDataSet();
+ gdcm::DataSet& dsInput = hreader.GetFile().GetDataSet();
gdcm::Attribute<0x54,0x11> series_number_att;
- series_number_att.SetFromDataSet(ds);
+ series_number_att.SetFromDataSet(dsInput);
int nbEnergyWindow = series_number_att.GetValue();
- std::cout << nbSlice << " " << nbEnergyWindow << std::endl;
+ const gdcm::DataElement &seqRotation = dsInput.GetDataElement(gdcm::Tag(0x54, 0x52));
+ gdcm::SmartPointer<gdcm::SequenceOfItems> sqiRotation = seqRotation.GetValueAsSQ();
+ gdcm::Item &itemRotation = sqiRotation->GetItem(1);
+ const gdcm::DataElement &deRotation = itemRotation.GetDataElement(gdcm::Tag(0x54, 0x53));
+ unsigned short SamplesPerPixelRotation = *((unsigned short *) deRotation.GetByteValue()->GetPointer());
+ std::stringstream s_SamplesPerPixelRotation;
+ s_SamplesPerPixelRotation << SamplesPerPixelRotation;
+ std::string p_StrRotation = s_SamplesPerPixelRotation.str();
+ int nbRotation = (int)atof(p_StrRotation.c_str());
int nbSlicePerEnergy = nbSlice / nbEnergyWindow;
if (nbSlicePerEnergy*nbEnergyWindow != nbSlice)
writer->SetInput( energyImage );
writer->SetImageIO( gdcmIO );
typename ReaderType::FileNamesContainer fileNamesOutput;
+ std::string extension = NumberToString(energy) + ".dcm";
+ std::string filename_out = m_ArgsInfo.output_arg + extension;
fileNamesOutput.push_back(filename_out);
writer->SetFileNames( fileNamesOutput );
writer->SetMetaDataDictionaryArray(&outputArray);
std::cerr << "Error: Exception thrown while writing the series!!" << std::endl;
std::cerr << excp << std::endl;
}
- }
-
-/*
-
- //Read the metadata informations in the mhd
- MetaObject tObj(3);
- tObj.AddUserField("NumberOfEnergyWindows", MET_STRING);
- tObj.AddUserField("NbOfHeads", MET_STRING);
- tObj.AddUserField("NbOfProjections", MET_STRING);
- tObj.Read(m_InputFileName.c_str());
- char *p_energyNumber, *p_headNumber, *p_rotationNumber;
- int energyNumber, headNumber, rotationNumber;
- p_energyNumber = static_cast<char*>(tObj.GetUserField("NumberOfEnergyWindows"));
- p_headNumber = static_cast<char*>(tObj.GetUserField("NbOfHeads"));
- p_rotationNumber = static_cast<char*>(tObj.GetUserField("NbOfProjections"));
- energyNumber = atoi(p_energyNumber);
- headNumber = atoi(p_headNumber);
- rotationNumber = atoi(p_rotationNumber);
- for (unsigned int i=0; i<energyNumber; ++i) {
- std::string tagEnergyName("EnergyWindow_");
- tagEnergyName += NumberToString(i).c_str();
- std::string tagEnergyThreshold(tagEnergyName), tagEnergyUphold(tagEnergyName);
- tagEnergyThreshold += "_threshold";
- tagEnergyUphold += "_uphold";
- tObj.AddUserField(tagEnergyName.c_str(), MET_STRING);
- tObj.AddUserField(tagEnergyThreshold.c_str(), MET_STRING);
- tObj.AddUserField(tagEnergyUphold.c_str(), MET_STRING);
- }
- tObj.Read(m_InputFileName.c_str());
- std::vector<char*> p_EnergyWindowName(energyNumber), p_EnergyWindowThreshold(energyNumber), p_EnergyWindowUphold(energyNumber);
- std::vector<int> energyWindowThreshold(energyNumber), energyWindowUphold(energyNumber);
- for (unsigned int i=0; i<energyNumber; ++i) {
- std::string tagEnergyName("EnergyWindow_");
- tagEnergyName += NumberToString(i).c_str();
- std::string tagEnergyThreshold(tagEnergyName), tagEnergyUphold(tagEnergyName);
- tagEnergyThreshold += "_threshold";
- tagEnergyUphold += "_uphold";
- p_EnergyWindowName[i] = static_cast<char*>(tObj.GetUserField(tagEnergyName.c_str()));
- p_EnergyWindowThreshold[i] = static_cast<char*>(tObj.GetUserField(tagEnergyThreshold.c_str()));
- p_EnergyWindowUphold[i] = static_cast<char*>(tObj.GetUserField(tagEnergyUphold.c_str()));
- energyWindowThreshold[i] = atoi(p_EnergyWindowThreshold[i]);
- energyWindowUphold[i] = atoi(p_EnergyWindowUphold[i]);
- }
-
-//TODO if the input size/spacing and dicom model ones are different
-
- // Create a new mhd image with the correct dicom order slices
- typename InputImageType::Pointer mhdCorrectOrder = InputImageType::New();
- mhdCorrectOrder->SetRegions(input->GetLargestPossibleRegion());
- mhdCorrectOrder->Allocate();
- unsigned int zAxis(0); //z value for the input mhd image
- for (unsigned int energy = 0; energy < energyNumber; ++energy) {
- for (unsigned int head = 0; head < headNumber; ++head) {
- for (unsigned int rotation = 0; rotation < rotationNumber; ++rotation) {
- std::cout << "Energy " << energy << " Head " << head << " Rotation " << rotation << std::endl;
-
- typename InputImageType::IndexType startIteratorIndexCorrectOrder; //pixel index of mhdCorrectOrder
- startIteratorIndexCorrectOrder[0] = 0;
- startIteratorIndexCorrectOrder[1] = 0;
- startIteratorIndexCorrectOrder[2] = rotation + head*rotationNumber + energy*headNumber;
-
- typename InputImageType::IndexType startIteratorIndexOriginalOrder; //pixel index of input mhd
- startIteratorIndexOriginalOrder[0] = 0;
- startIteratorIndexOriginalOrder[1] = 0;
- startIteratorIndexOriginalOrder[2] = head + energy*headNumber + rotation*energyNumber;
-
- typename InputImageType::SizeType regionSizeIterator;
- regionSizeIterator[0] = input->GetLargestPossibleRegion().GetSize()[0];
- regionSizeIterator[1] = input->GetLargestPossibleRegion().GetSize()[1];
- regionSizeIterator[2] = 1;
-
- typename InputImageType::RegionType regionIteratorCorrectOrder;
- regionIteratorCorrectOrder.SetSize(regionSizeIterator);
- regionIteratorCorrectOrder.SetIndex(startIteratorIndexCorrectOrder);
-
- typename InputImageType::RegionType regionIteratorOriginalOrder;
- regionIteratorOriginalOrder.SetSize(regionSizeIterator);
- regionIteratorOriginalOrder.SetIndex(startIteratorIndexOriginalOrder);
-
- itk::ImageRegionIterator<InputImageType> CorrectOrderIterator(mhdCorrectOrder,regionIteratorCorrectOrder);
- itk::ImageRegionIterator<InputImageType> OriginalOrderIterator(input,regionIteratorOriginalOrder);
- while(!CorrectOrderIterator.IsAtEnd()) {
- CorrectOrderIterator.Set(OriginalOrderIterator.Get());
- ++CorrectOrderIterator;
- ++OriginalOrderIterator;
- }
-
- ++zAxis;
- }
- }
- }
-
-
- // update output dicom keys/tags
- // string for distinguishing items inside sequence:
- const std::string ITEM_ENCAPSULATE_STRING("DICOM_ITEM_ENCAPSULATE");
- std::string tempString = ITEM_ENCAPSULATE_STRING + "01";
- typename ReaderSeriesType::DictionaryRawPointer inputDict = (*(readerSeries->GetMetaDataDictionaryArray()))[0];
- typename ReaderSeriesType::DictionaryArrayType outputArray;
- typename ReaderSeriesType::DictionaryRawPointer outputDict = new typename ReaderSeriesType::DictionaryType;
- CopyDictionary (*inputDict, *outputDict);
- itk::EncapsulateMetaData<std::string>(*outputDict, "0054|0011", NumberToString(energyNumber));
- itk::EncapsulateMetaData<std::string>(*outputDict, "0054|0021", NumberToString(headNumber));
- outputArray.push_back(outputDict);
-
- // Output directory and filenames
- typedef itk::ImageSeriesWriter<OutputImageType, OutputImageType> WriterSerieType;
- typename WriterSerieType::Pointer writerSerie = WriterSerieType::New();
- writerSerie->SetInput( mhdCorrectOrder );
- writerSerie->SetImageIO( gdcmIO );
- typename ReaderSeriesType::FileNamesContainer fileNamesOutput;
- fileNamesOutput.push_back(filename_out);
- writerSerie->SetFileNames( fileNamesOutput );
- writerSerie->SetMetaDataDictionaryArray(&outputArray);
-
-
- // Write
- try {
- if (m_ArgsInfo.verbose_flag)
- std::cout << writerSerie << std::endl;
- writerSerie->Update();
- } catch( itk::ExceptionObject & excp ) {
- std::cerr << "Error: Exception thrown while writing the series!!" << std::endl;
- std::cerr << excp << std::endl;
- }
-
-
- //Write sequence dicom tag with gdcm
- gdcm::Reader reader;
- reader.SetFileName( fileNamesOutput[0].c_str() );
- reader.Read();
- gdcm::File &file = reader.GetFile();
- gdcm::DataSet &ds2 = file.GetDataSet();
- const unsigned int ptr_len = 42;
- char *ptr = new char[ptr_len];
- memset(ptr,0,ptr_len);
-
- //Write rotation tag
- // Create a Sequence
- gdcm::SmartPointer<gdcm::SequenceOfItems> rotationSq = new gdcm::SequenceOfItems();
- rotationSq->SetLengthToUndefined();
- // Create a dataelement
- gdcm::DataElement rotationDE( gdcm::Tag(0x54, 0x53) );
- rotationDE.SetVR( gdcm::VR::US );
- char essai = (char)rotationNumber;
- char *p_essai = &essai;
- rotationDE.SetByteValue(p_essai, 1);
- // Create an item
- gdcm::Item rotationIt;
- rotationIt.SetVLToUndefined();
- gdcm::DataSet &rotationDS = rotationIt.GetNestedDataSet();
- rotationDS.Insert(rotationDE);
- rotationSq->AddItem(rotationIt);
- // Insert sequence into data set
- gdcm::DataElement rotationDEParent( gdcm::Tag(0x54, 0x52) );
- rotationDEParent.SetVR(gdcm::VR::SQ);
- rotationDEParent.SetValue(*rotationSq);
- rotationDEParent.SetVLToUndefined();
-
- //Write energy
- gdcm::DataElement energyDEParent( gdcm::Tag(0x54, 0x12) );
- energyDEParent.SetVR(gdcm::VR::SQ);
- // Create a Sequence
- gdcm::SmartPointer<gdcm::SequenceOfItems> energySq = new gdcm::SequenceOfItems();
- energySq->SetLengthToUndefined();
- for (unsigned int i=0; i<energyNumber; ++i) {
- gdcm::SmartPointer<gdcm::SequenceOfItems> energyThresholdSq = new gdcm::SequenceOfItems();
- energyThresholdSq->SetLengthToUndefined();
+ //Write sequence dicom tag with gdcm
+ gdcm::Reader reader;
+ reader.SetFileName( filename_out.c_str() );
+ reader.Read();
+ gdcm::File &fileOutput = reader.GetFile();
+ gdcm::DataSet &dsOutput = fileOutput.GetDataSet();
+ const unsigned int ptr_len = 42;
+ char *ptr = new char[ptr_len];
+ memset(ptr,0,ptr_len);
+ //Write rotation tag
+ // Create a Sequence
+ gdcm::SmartPointer<gdcm::SequenceOfItems> rotationSq = new gdcm::SequenceOfItems();
+ rotationSq->SetLengthToUndefined();
// Create a dataelement
- gdcm::DataElement energyThresholdDE( gdcm::Tag(0x54, 0x14) );
- gdcm::DataElement energyUpholdDE( gdcm::Tag(0x54, 0x15) );
- energyThresholdDE.SetVR( gdcm::VR::DS );
- energyUpholdDE.SetVR( gdcm::VR::DS );
- energyThresholdDE.SetByteValue(p_EnergyWindowThreshold[i], (uint32_t)strlen(p_EnergyWindowThreshold[i]));
- energyUpholdDE.SetByteValue(p_EnergyWindowUphold[i], (uint32_t)strlen(p_EnergyWindowUphold[i]));
+ gdcm::DataElement rotationDE( gdcm::Tag(0x54, 0x53) );
+ rotationDE.SetVR( gdcm::VR::US );
+ char essai = (char)nbRotation;
+ char *p_essai = &essai;
+ rotationDE.SetByteValue(p_essai, 1);
// Create an item
- gdcm::Item energyThresholdIt;
- energyThresholdIt.SetVLToUndefined();
- gdcm::DataSet &energyThresholdDS = energyThresholdIt.GetNestedDataSet();
- energyThresholdDS.Insert(energyThresholdDE);
- energyThresholdDS.Insert(energyUpholdDE);
- energyThresholdSq->AddItem(energyThresholdIt);
+ gdcm::Item rotationIt;
+ rotationIt.SetVLToUndefined();
+ gdcm::DataSet &rotationDS = rotationIt.GetNestedDataSet();
+ rotationDS.Insert(rotationDE);
+ rotationSq->AddItem(rotationIt);
// Insert sequence into data set
- gdcm::DataElement energyThresholdDEParent( gdcm::Tag(0x54, 0x13) );
- energyThresholdDEParent.SetVR(gdcm::VR::SQ);
- energyThresholdDEParent.SetValue(*energyThresholdSq);
- energyThresholdDEParent.SetVLToUndefined();
- // Create a dataelement
- gdcm::DataElement energyNameDE( gdcm::Tag(0x54, 0x18) );
- energyNameDE.SetVR( gdcm::VR::SH );
- energyNameDE.SetByteValue(p_EnergyWindowName[i], (uint32_t)strlen(p_EnergyWindowName[i]));
- // Create an item
- gdcm::Item energyIt;
- energyIt.SetVLToUndefined();
- gdcm::DataSet &energyDS = energyIt.GetNestedDataSet();
- energyDS.Insert(energyNameDE);
- energyDS.Insert(energyThresholdDEParent);
- energySq->AddItem(energyIt);
+ gdcm::DataElement rotationDEParent( gdcm::Tag(0x54, 0x52) );
+ rotationDEParent.SetVR(gdcm::VR::SQ);
+ rotationDEParent.SetValue(*rotationSq);
+ rotationDEParent.SetVLToUndefined();
+ //Write energy
+ gdcm::DataElement energyDEParent( gdcm::Tag(0x54, 0x12) );
+ energyDEParent.SetVR(gdcm::VR::SQ);
+ // Create a Sequence
+ gdcm::SmartPointer<gdcm::SequenceOfItems> energySq = new gdcm::SequenceOfItems();
+ energySq->SetLengthToUndefined();
+
+ //Read the caracteristics of this energy window : the name and the thresholds up and down
+ const gdcm::DataElement &seqWindow = dsInput.GetDataElement(gdcm::Tag(0x54, 0x12));
+ gdcm::SmartPointer<gdcm::SequenceOfItems> sqiWindow = seqWindow.GetValueAsSQ();
+ gdcm::Item &itemWindow = sqiWindow->GetItem(energy+1);
+ const gdcm::DataElement &deNameWindow = itemWindow.GetDataElement(gdcm::Tag(0x54, 0x18));
+ std::string nameWindow( deNameWindow.GetByteValue()->GetPointer(), deNameWindow.GetByteValue()->GetLength() );
+
+ const gdcm::DataElement &deThresholds = itemWindow.GetDataElement(gdcm::Tag(0x54, 0x13));
+ gdcm::SmartPointer<gdcm::SequenceOfItems> sqiThresholds = deThresholds.GetValueAsSQ();
+ for (unsigned int nbWindow=1; nbWindow <= sqiThresholds->GetNumberOfItems() ; ++nbWindow) {
+ gdcm::Item &itemThresholds = sqiThresholds->GetItem(nbWindow);
+ const gdcm::DataElement &deThresholdDown = itemThresholds.GetDataElement(gdcm::Tag(0x54, 0x14));
+ std::string p_StrThresholdDown( deThresholdDown.GetByteValue()->GetPointer(), deThresholdDown.GetByteValue()->GetLength() );
+ double thresholdDown = (double)atof(p_StrThresholdDown.c_str());
+ const gdcm::DataElement &deThresholdUp = itemThresholds.GetDataElement(gdcm::Tag(0x54, 0x15));
+ std::string p_StrThresholdUp( deThresholdUp.GetByteValue()->GetPointer(), deThresholdUp.GetByteValue()->GetLength() );
+ double thresholdUp = (double)atof(p_StrThresholdUp.c_str());
+
+ //Now write it !
+ gdcm::SmartPointer<gdcm::SequenceOfItems> energyThresholdSq = new gdcm::SequenceOfItems();
+ energyThresholdSq->SetLengthToUndefined();
+ // Create a dataelement
+ gdcm::DataElement energyThresholdDE( gdcm::Tag(0x54, 0x14) );
+ gdcm::DataElement energyUpholdDE( gdcm::Tag(0x54, 0x15) );
+ energyThresholdDE.SetVR( gdcm::VR::DS );
+ energyUpholdDE.SetVR( gdcm::VR::DS );
+ energyThresholdDE.SetByteValue(NumberToString(thresholdDown).c_str(), (uint32_t)strlen(NumberToString(thresholdDown).c_str()));
+ energyUpholdDE.SetByteValue(NumberToString(thresholdUp).c_str(), (uint32_t)strlen(NumberToString(thresholdUp).c_str()));
+ // Create an item
+ gdcm::Item energyThresholdIt;
+ energyThresholdIt.SetVLToUndefined();
+ gdcm::DataSet &energyThresholdDS = energyThresholdIt.GetNestedDataSet();
+ energyThresholdDS.Insert(energyThresholdDE);
+ energyThresholdDS.Insert(energyUpholdDE);
+ energyThresholdSq->AddItem(energyThresholdIt);
+ // Insert sequence into data set
+ gdcm::DataElement energyThresholdDEParent( gdcm::Tag(0x54, 0x13) );
+ energyThresholdDEParent.SetVR(gdcm::VR::SQ);
+ energyThresholdDEParent.SetValue(*energyThresholdSq);
+ energyThresholdDEParent.SetVLToUndefined();
+ // Create a dataelement
+ gdcm::DataElement energyNameDE( gdcm::Tag(0x54, 0x18) );
+ energyNameDE.SetVR( gdcm::VR::SH );
+ energyNameDE.SetByteValue(nameWindow.c_str(), (uint32_t)strlen(nameWindow.c_str()));
+ // Create an item
+ gdcm::Item energyIt;
+ energyIt.SetVLToUndefined();
+ gdcm::DataSet &energyDS = energyIt.GetNestedDataSet();
+ energyDS.Insert(energyNameDE);
+ energyDS.Insert(energyThresholdDEParent);
+ energySq->AddItem(energyIt);
+ }
+ // Insert sequence into data set
+ energyDEParent.SetValue(*energySq);
+ energyDEParent.SetVLToUndefined();
+ dsOutput.Insert(energyDEParent);
+ dsOutput.Insert(rotationDEParent);
+ gdcm::Writer w;
+ w.SetFile( fileOutput );
+ w.SetFileName( filename_out.c_str() );
+ w.Write();
}
- // Insert sequence into data set
- energyDEParent.SetValue(*energySq);
- energyDEParent.SetVLToUndefined();
- ds2.Insert(energyDEParent);
- ds2.Insert(rotationDEParent);
-
- gdcm::Writer w;
- w.SetFile( file );
- w.SetFileName( fileNamesOutput[0].c_str() );
- w.Write();
-*/
#else
std::cout << "Use GDCM2" << std::endl;
#endif