[Insight-users] translating image volume using affine transformation

michiel mentink michael.mentink at st-hughs.ox.ac.uk
Mon Nov 30 06:18:33 EST 2009


FIXED IT!

found the problem explained before (sorry, can't find that particular post
back, since I closed the window and don't
recall the search term.

The problem was that the direction of the image wasn't copied to the filter,
instead the filter got a identity direction.
I solved it by adding the lines:

 ImageType::DirectionType direction2 = reader->GetOutput()->GetDirection();
and subsituting direction by direction2 in
 filter->SetOutputDirection(direction2 );

Now it works!

greets, Michael

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

> 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/747af69b/attachment.htm>


More information about the Insight-users mailing list