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