[Insight-users] Unable to write dicom tag using GDCMImageIO
suicheng gu
gusuch at yahoo.com
Mon Oct 8 09:04:53 EDT 2012
Yes, I was trying to change the direction.
But I just used the resampledicom.cxx for input and output and modified the reconstruction part.
________________________________
From: Bill Lorensen <bill.lorensen at gmail.com>
To: suicheng gu <gusuch at yahoo.com>
Cc: ITK Users <insight-users at public.kitware.com>
Sent: Monday, October 8, 2012 12:44 AM
Subject: Re: [Insight-users] Unable to write dicom tag using GDCMImageIO
Please keep replies on the list...
Why are you resampling? I thought you wanted to change the directions.
On Sun, Oct 7, 2012 at 10:15 AM, suicheng gu <gusuch at yahoo.com> wrote:
attached, I just modified the ResampleDICOM.cxx file
>
>
>
>// Resample a DICOM study
>// Usage: ResampleDICOM InputDirectory OutputDirectory
>// xSpacing ySpacing zSpacing
>//
>// Example: ResampleDICOM CT CTResample 0 0 1.5
>// will read a series from the CT directory and create a
>// new series in the CTResample directory. The new series
>// will have the same x,y spacing as the input series,
but
>// will have a z-spacing of 1.5.
>//
>// Description:
>// ResampleDICOM resamples a DICOM series with user-specified
>// spacing. The program outputs a new DICOM series with a series
>// number set to 1001. All non-private DICOM tags are moved from the input
>// series to the output series. The Image Position Patient is adjusted
>// for each slice to reflect the z-spacing. The number of slices in
>// the output series may be larger or smaller due to changes in the
>// z-spacing. To retain the spacing for a given dimension, specify 0.
>//
>// The program progresses as follows:
>// 1) Read the input series
>// 2) Resample the series according to the user specified x-y-z
>// spacing.
>// 3) Create a MetaDataDictionary for each slice.
>// 4) Shift data to undo the effect of a rescale intercept by the
>// DICOM
reader
>// 5) Write the new DICOM series
>//
>
>#include "itkVersion.h"
>
>#include "itkImage.h"
>#include "itkMinimumMaximumImageFilter.h"
>
>#include "itkGDCMImageIO.h"
>#include "itkGDCMSeriesFileNames.h"
>#include "itkNumericSeriesFileNames.h"
>
>#include "itkImageSeriesReader.h"
>#include "itkImageSeriesWriter.h"
>
>#include "itkResampleImageFilter.h"
>#include "itkShiftScaleImageFilter.h"
>#include "itkMetaDataObject.h"
>#include "itkIdentityTransform.h"
>#include "itkLinearInterpolateImageFunction.h"
>#include "itkOrientImageFilter.h"
>#include <itksys/SystemTools.hxx>
>
>#if ITK_VERSION_MAJOR >= 4
>#include "gdcmUIDGenerator.h"
>#else
>#include "gdcm/src/gdcmFile.h"
>#include "gdcm/src/gdcmUtil.h"
>#endif
>
>#include <string>
>
>static void CopyDictionary (itk::MetaDataDictionary
&fromDict,
> itk::MetaDataDictionary &toDict);
>
>
>int main( int argc, char* argv[] )
>{
> // Validate input parameters
> if( argc < 2 )
> {
> std::cerr << "Usage: "
> << argv[0]
> << " InputDicomDirectory OutputDicomDirectory spacing_x spacing_y spacing_z"
> << std::endl;
> return EXIT_FAILURE;
> }
>
> const unsigned int InputDimension = 3;
> const unsigned int OutputDimension = 2;
>
> typedef signed short
PixelType;
>
> typedef itk::Image< PixelType, InputDimension >
> InputImageType;
> typedef itk::Image< PixelType, OutputDimension >
> OutputImageType;
> typedef itk::ImageSeriesReader< InputImageType >
> ReaderType;
> typedef itk::GDCMImageIO
> ImageIOType;
> typedef itk::GDCMSeriesFileNames
> InputNamesGeneratorType;
> typedef itk::NumericSeriesFileNames
> OutputNamesGeneratorType;
> typedef itk::IdentityTransform< double, InputDimension >
> TransformType;
> typedef itk::LinearInterpolateImageFunction< InputImageType, double >
> InterpolatorType;
> typedef itk::ResampleImageFilter< InputImageType, InputImageType >
> ResampleFilterType;
> typedef
itk::ShiftScaleImageFilter< InputImageType, InputImageType >
> ShiftScaleType;
> typedef itk::ImageSeriesWriter< InputImageType, OutputImageType >
> SeriesWriterType;
>
>////////////////////////////////////////////////
>// 1) Read the input series
>
> ImageIOType::Pointer gdcmIO = ImageIOType::New();
> InputNamesGeneratorType::Pointer inputNames = InputNamesGeneratorType::New();
> inputNames->SetInputDirectory( argv[1] );
>
> const ReaderType::FileNamesContainer & filenames =
> inputNames->GetInputFileNames();
>
> ReaderType::Pointer reader = ReaderType::New();
>
> reader->SetImageIO( gdcmIO );
>
reader->SetFileNames( filenames );
> try
> {
> reader->Update();
> }
> catch (itk::ExceptionObject &excp)
> {
> std::cerr << "Exception thrown while reading the series" << std::endl;
> std::cerr << excp << std::endl;
> return EXIT_FAILURE;
> }
>
>////////////////////////////////////////////////
>// 2) Resample the series
> InterpolatorType::Pointer interpolator = InterpolatorType::New();
>
> TransformType::Pointer transform = TransformType::New();
> transform->SetIdentity();
>
> const InputImageType::SpacingType& inputSpacing =
> reader->GetOutput()->GetSpacing();
> const InputImageType::RegionType& inputRegion =
>
reader->GetOutput()->GetLargestPossibleRegion();
> const InputImageType::SizeType& inputSize =
> inputRegion.GetSize();
> const double inputOrigin[3]={reader->GetOutput()->GetOrigin()[0],reader->GetOutput()->GetOrigin()[1],reader->GetOutput()->GetOrigin()[2]};
> // const double outputOrigin[3] = {inputOrigin[0],inputOrigin[1],inputOrigin[2]};
>
> std::cout << "The input series in directory " << argv[1]
> << " has " << filenames.size() << " files with spacing "
> << inputSpacing
> << std::endl;
>
> // Compute the size of the output. The user specifies a spacing on
> // the command line. If the spacing
is 0, the input spacing will be
> // used. The size (# of pixels) in the output is recomputed using
> // the ratio of the input and output sizes.
> InputImageType::SpacingType outputSpacing;
> outputSpacing[0] = inputSpacing[0];
> outputSpacing[1] = inputSpacing[2];
> outputSpacing[2] = inputSpacing[1];
>
> bool changeInSpacing = false;
> for (unsigned int i = 0; i < 3; i++)
> {
> if (outputSpacing[i] == 0.0)
> {
> outputSpacing[i] = inputSpacing[i];
> }
> else
> {
> changeInSpacing = true;
> }
> }
> InputImageType::SizeType outputSize;
> typedef InputImageType::SizeType::SizeValueType
SizeValueType;
> outputSize[0] = inputSize[0];
> outputSize[1] = inputSize[2];
> outputSize[2] = inputSize[1];
>
>
>
>////////////////////////////////////////////////
>// 3) Create a MetaDataDictionary for each slice.
>
> // Copy the dictionary from the first image and override slice
> // specific fields
> ReaderType::DictionaryRawPointer inputDict = (*(reader->GetMetaDataDictionaryArray()))[0];
> ReaderType::DictionaryArrayType outputArray;
>
> // To keep the new series in the same study as the original we need
> // to keep the same study UID. But we need new series and frame of
> // reference UID's.
>#if ITK_VERSION_MAJOR >= 4
> gdcm::UIDGenerator suid;
> std::string seriesUID = suid.Generate();
> gdcm::UIDGenerator fuid;
> std::string frameOfReferenceUID =
fuid.Generate();
>#else
> std::string seriesUID = gdcm::Util::CreateUniqueUID( gdcmIO->GetUIDPrefix());
> std::string frameOfReferenceUID = gdcm::Util::CreateUniqueUID( gdcmIO->GetUIDPrefix());
>#endif
> std::string studyUID;
> std::string sopClassUID;
> itk::ExposeMetaData<std::string>(*inputDict, "0020|000d", studyUID);
> itk::ExposeMetaData<std::string>(*inputDict, "0008|0016", sopClassUID);
> gdcmIO->KeepOriginalUIDOn();
>
> for (unsigned int f = 0; f < outputSize[2]; f++)
> {
> // Create a new dictionary for this slice
> ReaderType::DictionaryRawPointer dict = new ReaderType::DictionaryType;
>
> // Copy the dictionary from the first slice
> CopyDictionary (*inputDict, *dict);
>
> // Set the UID's for the study, series, SOP and frame of
reference
> itk::EncapsulateMetaData<std::string>(*dict,"0020|000d", studyUID);
> itk::EncapsulateMetaData<std::string>(*dict,"0020|000e", seriesUID);
> itk::EncapsulateMetaData<std::string>(*dict,"0020|0052", frameOfReferenceUID);
>
> #if ITK_VERSION_MAJOR >= 4
> gdcm::UIDGenerator sopuid;
> std::string sopInstanceUID = sopuid.Generate();
> #else
> std::string sopInstanceUID = gdcm::Util::CreateUniqueUID( gdcmIO->GetUIDPrefix());
> #endif
> itk::EncapsulateMetaData<std::string>(*dict,"0008|0018", sopInstanceUID);
> itk::EncapsulateMetaData<std::string>(*dict,"0002|0003", sopInstanceUID);
>
> // Change fields that are slice specific
> itksys_ios::ostringstream value;
> value.str("");
> value << f + 1;
>
> // Image Number
> itk::EncapsulateMetaData<std::string>(*dict,"0020|0013", value.str());
>
> // Series Description - Append new description to current series
> // description
> std::string oldSeriesDesc;
> itk::ExposeMetaData<std::string>(*inputDict, "0008|103e", oldSeriesDesc);
>
> value.str("");
> value << oldSeriesDesc
> << ": Resampled with pixel spacing "
> << outputSpacing[0] << ", "
> << outputSpacing[1] << ", "
> << outputSpacing[2];
> // This is an long string and there is a 64 character limit in the
> // standard
> unsigned lengthDesc = value.str().length();
>
> std::string seriesDesc( value.str(),
0,
> lengthDesc > 64 ? 64
> : lengthDesc);
> itk::EncapsulateMetaData<std::string>(*dict,"0008|103e", seriesDesc);
>
> // Series Number
> value.str("");
> value << 1001;
> itk::EncapsulateMetaData<std::string>(*dict,"0020|0011", value.str());
>
> // Derivation Description - How this image was derived
> value.str("");
> for (int i = 0; i < argc; i++)
> {
> value << argv[i] << " ";
> }
> value << ": " << ITK_SOURCE_VERSION;
>
> lengthDesc = value.str().length();
> std::string derivationDesc( value.str(), 0,
> lengthDesc > 1024 ?
1024
> : lengthDesc);
> itk::EncapsulateMetaData<std::string>(*dict,"0008|2111", derivationDesc);
>
> // Image orientation
> value.str("");
> value<<"1\\0\\0\\0\\0\\1";
> itk::EncapsulateMetaData<std::string>(*dict,"0020|0037", value.str());
> std::string orientation;
> itk::ExposeMetaData<std::string>(*dict, "0020|0037", orientation);
> //std::cout<<orientation<<std::endl;
> // Image Position Patient: This is calculated by computing the
> // physical coordinate of the first pixel in each slice.
> value.str("");
> value << inputOrigin[0] << "\\" <<(inputOrigin[1]+outputSpacing[2]*f)<<"\\"<< inputOrigin[2] ;
> //std::cout<<"position:
"<<value.str();
> itk::EncapsulateMetaData<std::string>(*dict,"0020|0032", value.str());
>
> // Slice Location: For now, we store the z component of the Image
> // Position Patient.
> value.str("");
> value << inputOrigin[1]+outputSpacing[2]*f;
> itk::EncapsulateMetaData<std::string>(*dict,"0020|1041", value.str());
>
> // Slice Thickness: For now, we store the z spacing
> value.str("");
> value << outputSpacing[2];
> itk::EncapsulateMetaData<std::string>(*dict,"0018|0050",
> value.str());
> // Spacing Between
Slices
> itk::EncapsulateMetaData<std::string>(*dict,"0018|0088",
> value.str());
>
> // Save the dictionary
> outputArray.push_back(dict);
> }
>
>////////////////////////////////////////////////
>// 4) Shift data to undo the effect of a rescale intercept by the
>// DICOM reader
> /*InputImageType::Pointer image = InputImageType::New();
> image->SetSpacing( outputSpacing );
> InputImageType::IndexType start;
> start[0] = 0;
> start[1] = 0;
> start[2] = 0;
> InputImageType::RegionType region;
> region.SetSize( outputSize );
> region.SetIndex( start );
> image->SetRegions( region );
> image->SetOrigin(inputOrigin);
>
> image->Allocate();
> //image->FillBuffer(0);
>
> PixelType *ppt = image->GetBufferPointer();
> for(int z = 0;z<outputSize[2];z++)
> {
> for(int y = 0;y<outputSize[1];y++)
> {
> for(int x = 0;x<outputSize[0];x++)
> {
> InputImageType::IndexType pixelIndex;
> pixelIndex[0] = x;
> pixelIndex[1] = z;
> pixelIndex[2] = y;
> PixelType pt = reader->GetOutput()->GetPixel(pixelIndex);
> *ppt = pt;
> ppt ++;
> }
> }
> }*/
> InputImageType::DirectionType direction;
> direction(0,0) = 1;
> direction(0,1) =
0;
> direction(0,2) = 0;
> direction(1,0) = 0;
> direction(1,1) = 1;
> direction(1,2) = 0;
> direction(2,0) = 0;
> direction(2,1) = 0;
> direction(2,2) = 1;
>
> //image->SetDirection(direction);
> //image->Update();
>
> /*
>
>itk::OrientImageFilter<InputImageType,InputImageType>::Pointer orienter =
>
>itk::OrientImageFilter<InputImageType,InputImageType>::New();
>
>orienter->UseImageDirectionOn();
>orienter->SetDesiredCoordinateOrientationToCoronal();
>orienter->SetInput(reader->GetOutput());
>
>orienter->Update();
>
*/
>itk::SpatialOrientation::ValidCoordinateOrientationFlags fileOrientation;
>
>//itk::ExposeMetaData<itk::SpatialOrientation::ValidCoordinateOrientationFlags>(reader->GetOutput()->GetMetaDataDictionary(),itk::ITK_CoordinateOrientation,fileOrientation);
>
>
>itk::OrientImageFilter<InputImageType,InputImageType>::Pointer orienter =
>
>itk::OrientImageFilter<InputImageType,InputImageType>::New();
>
>orienter->SetGivenCoordinateOrientation(itk::SpatialOrientation::ITK_COORDINATE_ORIENTATION_LPS); // deprecated
>
>orienter->SetDesiredCoordinateOrientation(itk::SpatialOrientation::ITK_COORDINATE_ORIENTATION_RIP);
>
>
>orienter->SetInput(reader->GetOutput());
>
>orienter->Update();
>
>////////////////////////////////////////////////
>// 5) Write the new DICOM series
>
> // Make the output directory and generate the file names.
> itksys::SystemTools::MakeDirectory( argv[2]
);
>
> // Generate the file names
> OutputNamesGeneratorType::Pointer outputNames = OutputNamesGeneratorType::New();
> std::string seriesFormat(argv[2]);
> seriesFormat = seriesFormat + "/" + "IM%d.dcm";
> outputNames->SetSeriesFormat (seriesFormat.c_str());
> outputNames->SetStartIndex (1);
> outputNames->SetEndIndex (outputSize[2]);
>
> //itk::MetaDataDictionary & dict1 = gdcmIO->GetMetaDataDictionary();
> //itk::EncapsulateMetaData<std::string>(dict1,"0020|0037", "1\\0\\0\\0\\0\\1");
>
> SeriesWriterType::Pointer seriesWriter = SeriesWriterType::New();
> seriesWriter->SetInput( orienter->GetOutput() );
> seriesWriter->SetImageIO( gdcmIO );
> seriesWriter->SetFileNames( outputNames->GetFileNames() );
>
seriesWriter->SetMetaDataDictionaryArray( &outputArray );
> try
> {
> seriesWriter->Update();
> }
> catch( itk::ExceptionObject & excp )
> {
> std::cerr << "Exception thrown while writing the series " << std::endl;
> std::cerr << excp << std::endl;
> return EXIT_FAILURE;
> }
> std::cout << "The output series in directory " << argv[2]
> << " has " << outputSize[2] << " files with spacing "
> << outputSpacing
> << std::endl;
> system("pause");
> return
EXIT_SUCCESS;
>}
>
>void CopyDictionary (itk::MetaDataDictionary &fromDict, itk::MetaDataDictionary &toDict)
>{
> typedef itk::MetaDataDictionary DictionaryType;
>
> DictionaryType::ConstIterator itr = fromDict.Begin();
> DictionaryType::ConstIterator end = fromDict.End();
> typedef itk::MetaDataObject< std::string > MetaDataStringType;
>
> while( itr != end )
> {
> itk::MetaDataObjectBase::Pointer entry = itr->second;
>
> MetaDataStringType::Pointer entryvalue =
> dynamic_cast<MetaDataStringType *>( entry.GetPointer() ) ;
> if( entryvalue )
> {
> std::string tagkey = itr->first;
> std::string tagvalue =
entryvalue->GetMetaDataObjectValue();
> itk::EncapsulateMetaData<std::string>(toDict, tagkey, tagvalue);
> }
> ++itr;
> }
>}
>
>
>From: Bill Lorensen <bill.lorensen at gmail.com>
>To: suicheng gu <gusuch at yahoo.com>
>Cc: "insight-users at itk.org" <insight-users at itk.org>
>Sent: Sunday, October 7, 2012 8:14 PM
>
>Subject: Re: [Insight-users] Unable to write dicom tag using GDCMImageIO
>
>
>
>Can you post the entire program please?
>
>
>On Sat, Oct 6, 2012 at 9:04 PM, suicheng gu <gusuch at yahoo.com> wrote:
>
>I used the following code to obtain the coronal view, but the orientation tag is still 1\0\0\0\1\0
>>And one more problem is that some dcm files are with size of 0.
>>
>>itk::OrientImageFilter<InputImageType,InputImageType>::Pointer orienter =
>>itk::OrientImageFilter<InputImageType,InputImageType>::New();
>>orienter->UseImageDirectionOn();
>>orienter->SetDesiredCoordinateOrientationToCoronal();
>>orienter->SetInput(reader->GetOutput());
>>orienter->Update();
>>
>>
>>From: Bill Lorensen <bill.lorensen at gmail.com>
>>To: suicheng gu <gusuch at yahoo.com>
>>Cc: "insight-users at itk.org" <insight-users at itk.org>
>>Sent: Sunday, October 7, 2012 4:33 AM
>>
>>Subject: Re: [Insight-users] Unable to write dicom tag using GDCMImageIO
>>
>>
>>You can use itkOrientImageFilter to change the direction of the volume. Then when you write the IDCOM with itk, the direction cosines will be correct.
>>
>>
>>On Sat, Oct 6, 2012 at 11:24 AM, suicheng gu <gusuch at yahoo.com> wrote:
>>
>>Thanks very much.
>>>what I want to do is to reconstruct the coronal and saggittal views and save them as dcm files.
>>>I can manually change the direction of the volume data
>>>If I can't change these tags, can I just create one dictionary? such as the example: ImageReadDicomSeriesWrite.cxx
>>>
>>>
>>>From: Bill Lorensen <bill.lorensen at gmail.com>
>>>To: suicheng gu <gusuch at yahoo.com>
>>>Cc: "insight-users at itk.org" <insight-users at itk.org>
>>>Sent: Saturday, October 6, 2012 11:07 PM
>>>Subject: Re: [Insight-users] Unable to write dicom tag using GDCMImageIO
>>>
>>>
>>>
>>>If I recall, tags that are related to orientation and spacing cannot be changed. You must change the orientation of the data before you write the dicom file. The writer will output the proper tag that goes with the data.
>>>
>>>
>>>On Sat, Oct 6, 2012 at 5:21 AM, suicheng gu <gusuch at yahoo.com> wrote:
>>>
>>>Hi,
>>>>
>>>>I was trying to chage the image orientation tag using
>>>>itk::EncapsulateMetaData<std::string>(dict1,"0020|0037", "1\\0\\0\\0\\0\\1");
>>>>But the image orientation was still 1\0\0\0\1\0
>>>>
>>>>
>>>>
>>>>I also tested the provided exmaple:
>>>>InsightToolkit-4.1.0\Examples\IO\DicomImageReadChangeHeaderWrite.cxx
>>>>
>>>>It works for some tags, such as the instance number: DicomImageReadChangeHeaderWrite $infile$ $outfile$ 0020|0013 10
>>>>
>>>>But it failed on the image orientation: DicomImageReadChangeHeaderWrite $infile$ $outfile$ 0020|0037 1\0\0\0\0\1
>>>>
>>>>
>>>>
>>>>Anyone can help?
>>>>Thanks in advance.Suicheng
>>>>
>>>>_____________________________________
>>>>
Powered by http://www.kitware.com/
>>>>
>>>>
>>>>Visit other Kitware open-source projects at
>>>>http://www.kitware.com/opensource/opensource.html
>>>>
>>>>Kitware offers ITK Training Courses, for more information visit:
>>>>http://www.kitware.com/products/protraining.php
>>>>
>>>>Please keep messages on-topic and check the ITK FAQ at:
>>>>http://www.itk.org/Wiki/ITK_FAQ
>>>>
>>>>Follow this link to subscribe/unsubscribe:
>>>>http://www.itk.org/mailman/listinfo/insight-users
>>>>
>>>>
>>>
>>>
>>>
>>>
>>>--
>>>Unpaid intern in BillsBasement at noware dot com
>>>
>>>
>>>
>>>
>>
>>
>>
>>--
>>Unpaid intern in BillsBasement at noware dot com
>>
>>
>>
>>
>
>
>
>--
>Unpaid intern in BillsBasement at noware dot com
>
>
>
>
--
Unpaid intern in BillsBasement at noware dot com
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20121008/7955a1c1/attachment-0001.htm>
More information about the Insight-users
mailing list