[Insight-users] problem with orientation of dicom output from segmentation
John Drozd
john.drozd at gmail.com
Fri Nov 27 12:42:31 EST 2009
Hi again,
I was at home and I just figured out how to provide you this link to my data
remotely,
so I provide it now instead of on the weekend:
Here is a link to my data (3D dicom volume):
input volume for my ConnectedThresholdFilter code:
http://www.apmaths.uwo.ca/~jdrozd/correctedsubject5.dcm
output volume from my ConnectedThresholdFilter code:
http://www.apmaths.uwo.ca/~jdrozd/outsubject5.dcm
If you need my CMakeLists.txt file, just let me know.
to run the code type:
./ConnectedThresholdImageFilter correctedsubject5.dcm outsubject5.dcm 103
142 95 17100 17300
and here is my code below:
Below is my uncommented code:
/*
to run type:
./ConnectedThresholdImageFilter correctedsubject5.dcm outsubject5.dcm 103
142 95 17100 17300
*/
#if defined(_MSC_VER)
#pragma warning ( disable : 4786 )
#endif
#ifdef __BORLANDC__
#define ITK_LEAN_AND_MEAN
#endif
#include "itkConnectedThresholdImageFilter.h"
#include "itkImage.h"
#include "itkCastImageFilter.h"
#include "itkCurvatureFlowImageFilter.h"
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"
#include "itkGDCMImageIO.h"
#include "itkVersion.h"
#include "itkOrientedImage.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 "itkIdentityTransform.h"
#include "itkLinearInterpolateImageFunction.h"
#include <itksys/SystemTools.hxx>
#include "gdcm/src/gdcmFile.h"
#include "gdcm/src/gdcmUtil.h"
#include <string>
int main( int argc, char *argv[])
{
if( argc < 7 )
{
std::cerr << "Missing Parameters " << std::endl;
std::cerr << "Usage: " << argv[0];
std::cerr << " inputImage outputImage seedX seedY seedZ lowerThreshold
upperThreshold" << std::endl;
return 1;
}
typedef float InternalPixelType;
const unsigned int Dimension = 3;
typedef itk::Image< InternalPixelType, Dimension > InternalImageType;
typedef signed short OutputPixelType;
typedef itk::Image< OutputPixelType, Dimension > OutputImageType;
typedef itk::Image< float, Dimension > OutputImageType2;
typedef itk::CastImageFilter< InternalImageType, OutputImageType >
CastingFilterType;
CastingFilterType::Pointer caster = CastingFilterType::New();
const unsigned int ImageDimension = 3;
typedef signed short PixelType;
typedef itk::Image< PixelType, ImageDimension > FixedImageType;
typedef itk::Image< float, ImageDimension > FloatImageType;
typedef itk::ImageFileReader< FixedImageType > ReaderType;
typedef itk::ImageFileWriter< OutputImageType > WriterType;
typedef itk::ImageFileWriter< FloatImageType > WriterType2;
ReaderType::Pointer reader = ReaderType::New();
WriterType::Pointer writer = WriterType::New();
WriterType2::Pointer writer2 = WriterType2::New();
typedef itk::GDCMImageIO ImageIOTypefixed;
ImageIOTypefixed::Pointer gdcmImageIOfixed = ImageIOTypefixed::New();
reader->SetImageIO( gdcmImageIOfixed );
typedef itk::GDCMImageIO ImageIOTypefixed2;
ImageIOTypefixed2::Pointer gdcmImageIOfixed2 = ImageIOTypefixed2::New();
reader->SetFileName( argv[1] );
reader->Update();
typedef itk::CurvatureFlowImageFilter< InternalImageType,
InternalImageType >
CurvatureFlowImageFilterType;
CurvatureFlowImageFilterType::Pointer smoothing =
CurvatureFlowImageFilterType::New();
typedef itk::ConnectedThresholdImageFilter< InternalImageType,
InternalImageType > ConnectedFilterType;
ConnectedFilterType::Pointer connectedThreshold =
ConnectedFilterType::New();
typedef signed short InputAPixelType;
typedef float OutputBPixelType;
typedef itk::Image< InputAPixelType, 3 > InputAImageType;
typedef itk::Image< OutputBPixelType, 3 > OutputBImageType;
typedef itk::CastImageFilter< InputAImageType, OutputBImageType >
CastFilterType;
CastFilterType::Pointer castFilter = CastFilterType::New();
castFilter->SetInput( reader->GetOutput() );
connectedThreshold->SetInput( castFilter->GetOutput() );
caster->SetInput( connectedThreshold->GetOutput() );
smoothing->SetNumberOfIterations( 20 ); //was 5
smoothing->SetTimeStep( 0.125 );
const InternalPixelType lowerThreshold = atof( argv[6] );
const InternalPixelType upperThreshold = atof( argv[7] );
connectedThreshold->SetLower( lowerThreshold );
connectedThreshold->SetUpper( upperThreshold );
connectedThreshold->SetReplaceValue( 255 );
InternalImageType::IndexType index;
index[0] = atoi( argv[3] );
index[1] = atoi( argv[4] );
//added
index[2] = atoi( argv[5] );
std::cout << index << std::endl;
// Software Guide : BeginCodeSnippet
connectedThreshold->SetSeed( index );
//obtain a 5 x 5 bounding region of seeds
int ii, jj, kk;
ii = index[0];
jj = index[1];
kk = index[2];
for (int i = ii; i < ii + 5; i++)
for (int j = jj; j < jj + 5; j++)
for (int k = kk; k < kk + 5; k++)
{
index[0] = i;
index[1] = j;
index[2] = k;
connectedThreshold->AddSeed( index );
}
for (int i = ii; i > ii - 5; i--)
for (int j = jj; j > jj - 5; j--)
for (int k = kk; k > kk - 5; k--)
{
index[0] = i;
index[1] = j;
index[2] = k;
connectedThreshold->AddSeed( index );
}
connectedThreshold->Print(std::cout,17100);
typedef itk::MetaDataDictionary DictionaryType;
DictionaryType inputdict = reader->GetMetaDataDictionary();
writer->SetMetaDataDictionary( inputdict );
writer->SetFileName( argv[2] );
writer->SetInput( caster->GetOutput() );
try
{
writer->Update();
}
catch( itk::ExceptionObject & excep )
{
std::cerr << "Exception caught !" << std::endl;
std::cerr << excep << std::endl;
}
return 0;
On Fri, Nov 27, 2009 at 12:31 PM, John Drozd <john.drozd at gmail.com> wrote:
> Hi all,
>
> I originally have a dicom series of 2d slices which I had converted to a 3D
> dicom volume for processing.
> I'll try processing the dicom series and then write the output to a dicom
> series of 2d slices, to see if I can circumvent this problem.
> I'll let you know if this works?
>
> I'll provide a link to my data as well, but it probably won't be until the
> weekend.
>
> Thanks,
> john
>
> On Fri, Nov 27, 2009 at 12:23 PM, Mathieu Malaterre <
> mathieu.malaterre at gmail.com> wrote:
>
>> This is a post ITK 3.16 feature I am afraid:
>>
>> http://www.cmake.org/Bug/view.php?id=7748
>>
>> (0017626)
>> Mathieu.Malaterre (developer)
>> 2009-09-19 06:51
>>
>> This will be post ITK 3.16.
>>
>> Thanks for the patch ! Sorry for being so long to fix it.
>>
>> On Fri, Nov 27, 2009 at 6:11 PM, Bill Lorensen <bill.lorensen at gmail.com>
>> wrote:
>> > John,
>> >
>> > 3.16 will not help. This looks like a bug for sure.
>> >
>> > Can you provide a link to a dataset that we can use to reproduce the
>> > problem? We don't have much experience with 3D dicom files I'm afraid.
>> >
>> > Bill
>> >
>> > On Fri, Nov 27, 2009 at 11:59 AM, John Drozd <john.drozd at gmail.com>
>> wrote:
>> >> Hi Bill,
>> >>
>> >> I am currently using the ITK 3.14 that is in the 3D Slicer 3.4
>> >> directory Slicer3-lib/Insight
>> >>
>> >> If I upgrade to ITK 3.16, will this solve my problem?
>> >>
>> >> I downloaded ITK 3.16 and noticed using the diff command that
>> >> itkGDCMIO.h has been modified from what it was in ITK 3.14.
>> >>
>> >> But when I ran diff on gdcmOrientation.cxx in both versions of ITK,
>> >> diff showed no difference.
>> >>
>> >> John
>> >>
>> >> On Friday, November 27, 2009, Bill Lorensen <bill.lorensen at gmail.com>
>> wrote:
>> >>> John,
>> >>>
>> >>> What version of itk are you using? Are you using the gdcm 1.x that is
>> >>> in the itk source tree or an external gdcm 2.x version?
>> >>>
>> >>> Bill
>> >>>
>> >>> On Thu, Nov 26, 2009 at 3:22 PM, John Drozd <john.drozd at gmail.com>
>> wrote:
>> >>>> Hi Luis,
>> >>>>
>> >>>> I was getting a segmentation fault because I had forgotten to change
>> the
>> >>>> Dimension from 2 to 3 in DicomReadPrintTag.cxx.
>> >>>>
>> >>>> For the input image, I have:
>> >>>> (0020|0037) Image Orientation (Patient) = 0.0\0.0\-1.0\0.0\1.0\0.0
>> >>>>
>> >>>> For the output image, I have:
>> >>>> (0020|0037) Image Orientation (Patient) =
>> >>>> 0.000000\0.000000\1.000000\0.000000\1.000000\0.000000
>> >>>>
>> >>>> There is a difference. Since the program is changing the orientation
>> when
>> >>>> processed through the pipeline, is the best approach to manually
>> change the
>> >>>> orientation tag of the dictionary and then copy the revised
>> dictionary to
>> >>>> the output file?
>> >>>>
>> >>>> john
>> >>>>
>> >>>> Below are the full tags:
>> >>>>
>> >>>> [jdrozd at trumpet DicomImageReadPrintTags]$ ./DicomImageReadPrintTags
>> >>>> correctedsubject5.dcm
>> >>>> (0002|0000) Group Length =
>> >>>> 214
>> >>>> (0002|0001) File Meta Information Version =
>> >>>> AAE=
>> >>>> (0002|0002) Media Storage SOP Class UID =
>> >>>> 1.2.840.10008.5.1.4.1.1.2
>> >>>> (0002|0003) Media Storage SOP Instance UID =
>> >>>> 1.2.826.0.1.3680043.2.1125.1.10607669833050788267094246636093811
>> >>>> (0002|0010) Transfer Syntax UID =
>> >>>> 1.2.840.10008.1.2.1
>> >>>> (0002|0012) Implementation Class UID =
>> >>>> 147.144.143.155
>> >>>> (0002|0013) Implementation Version Name = ITK/GDCM
>> >>>> 1.2.4
>> >>>> (0002|0016) Source Application Entity Title =
>> >>>> NOTSPECIFIED
>> >>>> (0008|0008) Image Type =
>> >>>> DERIVED\PRIMARY
>> >>>> (0008|0012) Instance Creation Date =
>> >>>> 20091029
>> >>>> (0008|0013) Instance Creation Time =
>> >>>> 135224
>> >>>> (0008|0016) SOP Class UID =
>> >>>> 1.2.840.10008.5.1.4.1.1.2
>> >>>> (0008|0018) SOP Instance UID =
>> >>>> 1.2.826.0.1.3680043.2.1125.1.10607669833050788267094246636093811
>> >>>> (0008|0020) Study Date =
>> >>>> 20081030
>> >>>> (0008|0030) Study Time =
>> >>>> 164348.940
>> >>>> (0008|0050) Accession Number
>> >>>> =
>> >>>> (0008|0060) Modality =
>> >>>> CT
>> >>>> (0008|0064) Conversion Type =
>> >>>> WSD
>> >>>> (0008|0070) Manufacturer =
>> >>>> Manifacturer
>> >>>> (0008|0080) Institution Name = GDCM
>> >>>> Hospital
>> >>>> (0008|0090) Referring Physician's Name = Refering
>> >>>> Phisician
>> >>>> (0010|0010) Patient's Name = Patient
>> >>>> 188858520
>> >>>> (0010|0020) Patient ID =
>> >>>> 1747233212
>> >>>> (0010|0030) Patient's Birth Date =
>> >>>> 19500101
>> >>>> (0010|0040) Patient's Sex =
>> >>>> M
>> >>>> (0018|0088) Spacing Between Slices =
>> >>>> 1.207500
>> >>>> (0018|1164) Imager Pixel Spacing =
>> >>>> 0.945750\0.945750
>> >>>> (0020|000d) Study Instance UID =
>> >>>> 1.2.826.0.1.3680043.2.1125.1.53653479342656887425609263146205769
>> >>>> (0020|000e) Series Instance UID =
>> >>>> 1.2.826.0.1.3680043.2.1125.1.47155241092896696198844159625679986
>> >>>> (0020|0010) Study ID =
>> >>>> 1533117581
>> >>>> (0020|0011) Series Number =
>> >>>> 2135500125
>> >>>> (0020|0013) Instance Number =
>> >>>> 165
>> >>>> (0020|0020) Patient Orientation =
>> >>>> L\P
>> >>>> (0020|0032) Image Position (Patient) =
>> >>>> 199.237496852874\0.0\0.0
>> >>>> (0020|0037) Image Orientation (Patient) =
>> >>>> 0.0\0.0\-1.0\0.0\1.0\0.0
>> >>>> (0020|4000) Image Comments = NOT FOR CLINICAL
>> >>>> USE
>> >>>> (0028|0002) Samples per Pixel =
>> >>>> 1
>> >>>> (0028|0004) Photometric Interpretation =
>> >>>> MONOCHROME2
>> >>>> (0028|0008) Number of Frames =
>> >>>> 166
>> >>>> (0028|0010) Rows =
>> >>>> 256
>> >>>> (0028|0011) Columns =
>> >>>> 256
>> >>>> (0028|0030) Pixel Spacing =
>> >>>> 0.945750\0.945750
>> >>>> (0028|0034) Pixel Aspect Ratio =
>> >>>> 1\1
>> >>>> (0028|0100) Bits Allocated =
>> >>>> 16
>> >>>> (0028|0101) Bits Stored =
>> >>>> 16
>> >>>> (0028|0102) High Bit =
>> >>>> 15
>> >>>> (0028|0103) Pixel Representation =
>> >>>> 0
>> >>>> (0028|1052) Rescale Intercept =
>> >>>> 0.0
>> >>>> (0028|1053) Rescale Slope =
>> >>>> 1.0
>> >>>> (0028|1054) Rescale Type =
>> >>>> US
>> >>>> (7fe0|0000) Group Length =
>> >>>> 21757964
>> >>>> Patient's Name (0010|0010) is: Patient
>> >>>> 188858520
>> >>>> Performing Physician's Name (0008|1050): (No Value Found in
>> >>>> File)
>> >>>> PixelType:
>> >>>> scalar
>> >>>> Component Type: unsigned_short
>> >>>>
>> >>>> [jdrozd at trumpet DicomImageReadPrintTags]$ ./DicomImageReadPrintTags
>> >>>> outsubject5.dcm
>> >>>> (0002|0000) Group Length =
>> >>>> 194
>> >>>> (0002|00
>> >>
>> > _____________________________________
>> > Powered by 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.html
>> >
>> > 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
>> >
>>
>>
>>
>> --
>> Mathieu
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20091127/e05fe929/attachment-0001.htm>
More information about the Insight-users
mailing list