[Insight-users] translating image volume using affine transformation

michiel mentink michael.mentink at st-hughs.ox.ac.uk
Mon Nov 30 05:10:49 EST 2009


To check this, I also followed Luiz Ibanez' suggestion to John Drozd in the
message "problem with orientation of
dicom output from segmentation.

I pasted in the following lines into the codes:

std::cout << reader->GetOutput()->GetDirection() << std::endl;
std::cout << filter->GetOutput()->GetDirection() << std::endl;

and got the following result:

0.997256 -0.0188216 -0.0716039
0.0740362 0.253523 0.964492
-8.34627e-10 -0.967146 0.254221

1 0 0
0 1 0
0 0 1

Does this mean for certain that the filter changes the axis' orientation?

Kind regards,

Michael


On Mon, Nov 30, 2009 at 9:53 AM, michiel mentink <
michael.mentink at st-hughs.ox.ac.uk> wrote:

> might this be the same problem as in
> problem with orientation of dicom output from segmentation where the dicom
> 3d volume is changed in orientation?
>
> Kind regards, Michael
>
>
> On Mon, Nov 30, 2009 at 9:42 AM, michiel mentink <
> michael.mentink at st-hughs.ox.ac.uk> wrote:
>
>> Dear Bill,
>>
>> I indeed didn't define the origin correctly. However, now there seems to
>> be a more
>> serious problem.
>>
>> Since I'm not sure I can upload screenshots to this mailing list, I have
>> put a
>> ITKsnap screenshot of the original dcm volume view, together with output
>> after
>> the affine filter and code on this site:
>>
>>
>> http://sites.google.com/site/michielmentink/programming/itk/image-translate
>>
>> You can see there that before and after the affine filter, using identity
>> transformation,
>> the filesize, number of pixels, and pixel sizes are identical.
>>
>> There is a small bit of the original volume visible in the resulting dicom
>> volume, but
>> it contains only image data on a eight slices and the slices are deformed.
>>
>>
>> What I can make of the data, it seems that the filter has rotated the
>> dicom volume
>> so that the 'coronal' (seeing)from the front) view now shows the 'axial'
>> view (seeing
>> from above).
>> Also, the image is flipped in both X and Y axis.
>>
>> If I change the Z origin to -69 instead of 69, most image slices actually
>> contain the
>> image data, but still the views have been changed completely, like
>> mentioned before.
>>
>>
>> What is going on?
>> What I really want to do is load in a dicom volume, translate it in the X,
>> Y, Z direction
>> and save it. Normally, the affine filter should be used for this, right?
>>
>> For completeness, I've pasted the code again.
>>
>> Kind regards, Michael
>>
>>
>>
>>
>>
>> On Fri, Nov 27, 2009 at 4:05 PM, Bill Lorensen <bill.lorensen at gmail.com>wrote:
>>
>>> What is the origin, spacing and direction of the input volume?
>>>
>>> On Fri, Nov 27, 2009 at 9:59 AM, michiel mentink
>>> <michael.mentink at st-hughs.ox.ac.uk> wrote:
>>> > Dear all,
>>> >
>>> > I'm trying to make a simple function that translates a dicom volume by
>>> a
>>> > certain amount.
>>> > It first opens all dicom slides of an mri scan.
>>> >
>>> > To do this, I pasted DicomSeriesReadImageWrite2.cxx and
>>> > ResampleImageFilter.cxx and adjusted
>>> > the latter so it accepts 3 dimensions instead of two and uses integer
>>> > instead of unsigned char.
>>> >
>>> > The function compiles and a .vtk volume is written. However, there's no
>>> > image data in the file.
>>> > It does display all pixels with value 100, as I set by:
>>> > filter->SetDefaultPixelValue( 100 );
>>> >
>>> > questions:
>>> > 1) am I too naive to think that the ResampleImageFilter can easily be
>>> > extended from 2 to 3 dimensions?
>>> > 2) if not, what did I do wrong?
>>> >
>>>
>>
>>
>> #if defined(_MSC_VER)
>> #pragma warning ( disable : 4786 )
>> #endif
>>
>> #ifdef __BORLANDC__
>> #define ITK_LEAN_AND_MEAN
>> #endif
>>
>> #include "itkOrientedImage.h"
>> #include "itkGDCMImageIO.h"
>> #include "itkGDCMSeriesFileNames.h"
>> #include "itkImageSeriesReader.h"
>> #include "itkImageFileWriter.h"
>> #include "itkResampleImageFilter.h"
>>
>> #include "itkAffineTransform.h"
>> #include "itkNearestNeighborInterpolateImageFunction.h"
>> // Software Guide : EndCodeSnippet
>>
>> int main( int argc, char* argv[] )
>> {
>>
>>   if( argc < 3 )
>>     {
>>     std::cerr << "Usage: " << std::endl;
>>     std::cerr << argv[0] << " DicomDirectory  outputFileName
>> [seriesName]"
>>               << std::endl;
>>     return EXIT_FAILURE;
>>     }
>>
>>   typedef signed short    PixelType;
>>   const unsigned int      Dimension = 3;
>>
>>   typedef itk::OrientedImage< PixelType, Dimension >         ImageType;
>>
>>   typedef itk::ImageSeriesReader< ImageType >        ReaderType;
>>   ReaderType::Pointer reader = ReaderType::New();
>>
>>   typedef itk::GDCMImageIO       ImageIOType;
>>   ImageIOType::Pointer dicomIO = ImageIOType::New();
>>
>>   reader->SetImageIO( dicomIO );
>>
>>   typedef itk::GDCMSeriesFileNames NamesGeneratorType;
>>   NamesGeneratorType::Pointer nameGenerator = NamesGeneratorType::New();
>>
>>   nameGenerator->SetUseSeriesDetails( true );
>>   nameGenerator->AddSeriesRestriction("0008|0021" );
>>
>>   nameGenerator->SetDirectory( argv[1] );
>>
>>   try
>>     {
>>     std::cout << std::endl << "The directory: " << std::endl;
>>     std::cout << std::endl << argv[1] << std::endl << std::endl;
>>     std::cout << "Contains the following DICOM Series: ";
>>     std::cout << std::endl << std::endl;
>>
>>     typedef std::vector< std::string >    SeriesIdContainer;
>>
>>     const SeriesIdContainer & seriesUID = nameGenerator->GetSeriesUIDs();
>>
>>     SeriesIdContainer::const_iterator seriesItr = seriesUID.begin();
>>     SeriesIdContainer::const_iterator seriesEnd = seriesUID.end();
>>     while( seriesItr != seriesEnd )
>>       {
>>       std::cout << seriesItr->c_str() << std::endl;
>>       seriesItr++;
>>       }
>>
>>     std::string seriesIdentifier;
>>
>>     if( argc > 3 ) // If no optional series identifier
>>       {
>>       seriesIdentifier = argv[3];
>>       }
>>     else
>>       {
>>       seriesIdentifier = seriesUID.begin()->c_str();
>>       }
>>
>>     std::cout << std::endl << std::endl;
>>     std::cout << "Now reading series: " << std::endl << std::endl;
>>     std::cout << seriesIdentifier << std::endl;
>>     std::cout << std::endl << std::endl;
>>
>>     typedef std::vector< std::string >   FileNamesContainer;
>>     FileNamesContainer fileNames;
>>
>>     fileNames = nameGenerator->GetFileNames( seriesIdentifier );
>>     reader->SetFileNames( fileNames );
>>     try
>>       {
>>       reader->Update();
>>       }
>>     catch (itk::ExceptionObject &ex)
>>       {
>>       std::cout << ex << std::endl;
>>       return EXIT_FAILURE;
>>       }
>>     }
>>   catch (itk::ExceptionObject &ex)
>>     {
>>     std::cout << ex << std::endl;
>>     return EXIT_FAILURE;
>>     }
>> //  return EXIT_SUCCESS;
>>
>>
>> //////////////////////////////////////////////////////////////////////////////////////
>> // FROM HERE NEW CODE
>>
>>
>> /////////////////////////////////////////////////////////////////////////////////////
>>
>> //  const     unsigned int   Dimension = 3;
>> //  typedef   unsigned char  InputPixelType;
>> //  typedef   unsigned char  OutputPixelType;
>>   typedef itk::Image< PixelType, Dimension >   InputImageType;
>>    typedef itk::Image< PixelType, Dimension >   OutputImageType;
>>   typedef itk::ImageFileWriter< ImageType >  WriterType;
>>   WriterType::Pointer writer = WriterType::New();
>>
>>   writer->SetFileName( argv[2] );
>>
>>   typedef itk::ResampleImageFilter<ImageType,ImageType> FilterType;
>>   FilterType::Pointer filter = FilterType::New();
>>
>>   typedef itk::AffineTransform< double, Dimension >  TransformType;
>>   TransformType::Pointer transform = TransformType::New();
>>   filter->SetTransform( transform );
>>
>>
>>   typedef itk::NearestNeighborInterpolateImageFunction<ImageType, double
>> >  InterpolatorType;
>>   InterpolatorType::Pointer interpolator = InterpolatorType::New();
>>   filter->SetInterpolator( interpolator );
>>
>>   filter->SetDefaultPixelValue( 100 ); // if pixel in output image not
>> defined, what value?
>>
>>   double spacing[ Dimension ];
>>   spacing[0] = 0.29; // pixel spacing in millimeters along X
>>   spacing[1] = 0.29; // pixel spacing in millimeters along Y
>>   spacing[2] = 3.3;  // pixel spacing in millimeters along Z
>>
>>
>>   filter->SetOutputSpacing( spacing );
>>
>>   double origin[ Dimension ];
>>   origin[0] = -174.22;  // X space coordinate of origin
>>   origin[1] = -125.00;  // Y space coordinate of origin
>>   origin[2] =   69.37;  // Z space coordinate of origin
>>
>>
>>   filter->SetOutputOrigin( origin );
>>   ImageType::DirectionType direction;
>>   direction.SetIdentity();
>>   filter->SetOutputDirection( direction );
>>   ImageType::SizeType   size;
>>
>>   size[0] = 560;  // number of pixels along X
>>   size[1] = 560;  // number of pixels along Y
>>   size[2] = 26;
>>
>>
>>   filter->SetSize( size );
>>   filter->SetInput( reader->GetOutput() );
>>
>>   writer->SetInput( filter->GetOutput() );
>>   writer->Update();
>>
>>   std::cout << "Writing the image" << std::endl << std::endl;
>>
>>
>>   /* (this code not executed)
>>
>>   TransformType::OutputVectorType translation;
>>   translation[0] = 0;  // X translation in millimeters
>>   translation[1] = 0;  // Y translation in millimeters
>>   translation[2] = 0;  // Z translation in millimeters
>>
>>
>>   transform->Translate( translation );
>>   // Software Guide : EndCodeSnippet
>>
>>    writer->SetInput( filter->GetOutput() );
>>
>>     std::cout  << "Writing the image as " << std::endl << std::endl;
>>     std::cout  << argv[2] << std::endl << std::endl;
>>     writer->Update();
>>     */
>>
>> }
>>
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20091130/2a533cd6/attachment.htm>


More information about the Insight-users mailing list