(0014233)
Mathieu Malaterre (developer)
2008-11-30 08:56

Please attach demo, thanks.

(0014264)
Andreas Schuh (reporter)
2008-12-01 18:29

The Demo is quite simple. Try it with a volumetric DICOM file with the DICOM entry (0020,0037) Image Orientation (Patient) be 0/0/-1/1/0/0 for instance.

GDCMImageIO::InternalReadImageInformation() reads in the patient orientation correctly and sets it with ImageIOBase::SetDirection() for all three axes. The ImageFileReader then inverts this direction matrix by just setting the row vectors of ImageIOBase::m_Direction as columns of the image direction (which is the inverse patient orientation, the mapping from image space to physical space). Accordingly, ImageFileWriter takes the direction matrix of the image (ImageBase::GetDirection()) and sets the columns of this matrix as rows of the m_Direction member of the ImageIOBase subclass (in our case GDCMImageIO) via ImageIOBase::SetDirection(). Till that point everything is fine. But then in GDCMImageIO::Write() the following code writes the wrong elements of m_Direction as values of the (0020,0037) header entry:

if( !header->GetValEntry(0x0020,0x0037 ) )
{
  str.str("");
  str << m_Direction[0][0] << "\\"
      << m_Direction[1][0] << "\\";
   /*
    * This is where the 3rd component of the direction is being lost
    * ITK mechanism does not support 2D image, placed in 3D world...
    */
  if( m_Direction.size() == 3 )
  {
    str << m_Direction[2][0] << "\\";
  }
  else
  {
    str << 0. << "\\";
  }
  str << m_Direction[0][1] << "\\"
      << m_Direction[1][1] << "\\";
  if( m_Direction.size() == 3 )
  {
    str << m_Direction[2][1];
  }
  else
  {
    str << 0.;
  }
  header->InsertValEntry(str.str(),0x0020,0x0037);
}

Actually, (m_Direction[0][0], m_Direction[0][1], m_Direction[0][2]) is the direction of the first axis and (m_Direction[1][0], m_Direction[1][1], m_Direction[1][2]) the one of the second axis. Exactly how ImageFileWriter sets the directions via ImageIOBase::SetDirection().

(0014265)
Andreas Schuh (reporter)
2008-12-01 18:38

To resolve some confusion:

Patient Orientation: 0/0/-1/1/0/0

GDCMImageIO::m_Direction:

0 0 -1
1 0 0
0 -1 0

Image::m_Direction:

0 1 0
0 0 -1
-1 0 0

The conversion/inversion from one to another is done by ImageFileReader, ImageFileWriter.

(0014889)
Andreas Schuh (reporter)
2009-02-10 18:58

Hi Mathieu,

maybe you still don't believe me :-) ... Let me explain it another way ...

itk::GDCMImageIO::InternalReadImageInformation() sets the m_Direction member of the itk::ImageIOBase class doing the following:

  m_Direction.resize(3);
  this->SetDirection(0, rowDirection);
  this->SetDirection(1, columnDirection);
  this->SetDirection(2, sliceDirection);

So far so good. Having a look at itk::ImageIOBase::SetDirection() reveals the following:

  m_Direction[0] stores the rowDirection
  m_Direction[1] stores the columnDirection
  m_Direction[2] stores the sliceDirection

So,

  rowDirection = [m_Direction[0][0], m_Direction[0][1], m_Direction[0][2]]
  columnDirection = [m_Direction[1][0], m_Direction[1][1], m_Direction[1][2]]
  sliceDirection = [m_Direction[2][0], m_Direction[2][1], m_Direction[2][2]]

Any objections? Don't think so. So proceeding to the BUGGY code in itk::GDCMImageIO::Write().

/////////////////////////////////
// ATTENTION: THIS IS WRONG!!! //
/////////////////////////////////

  // Handle Direction = Image Orientation Patient
  // Same comment as above, if user tell us what the Orientation is, we should not try
  // to set if from the Image as we might have lost some information
  if( !header->GetValEntry(0x0020,0x0037 ) )
    {
    str.str("");
    str << m_Direction[0][0] << "\\"
      << m_Direction[1][0] << "\\";
    /*
     * This is where the 3rd component of the direction is being lost
     * ITK mechanism does not support 2D image, placed in 3D world...
     */
    if( m_Direction.size() == 3 )
      {
      str << m_Direction[2][0] << "\\";
      }
    else
      {
      str << 0. << "\\";
      }
    str << m_Direction[0][1] << "\\"
      << m_Direction[1][1] << "\\";
    if( m_Direction.size() == 3 )
      {
      str << m_Direction[2][1];
      }
    else
      {
      str << 0.;
      }
    header->InsertValEntry(str.str(),0x0020,0x0037); // Image Orientation (Patient)
    }

As you see, [m_Direction[0][0], m_Direction[1][0], m_Direction[2][0]] is written as row direction and [m_Direction[0][1], m_Direction[1][1], m_Direction[2][1]] as column direction ... but this is NOT right as we just saw by looking at itk::GDCMImageIO::InternalReadImageInformation().

Instead the code in itk::GDCMImageIO::Write() should look like this:

///////////////////////
// HERE'S THE BUGFIX //
///////////////////////

  // Handle Direction = Image Orientation Patient
  // Same comment as above, if user tell us what the Orientation is, we should not try
  // to set if from the Image as we might have lost some information
  if( !header->GetValEntry(0x0020,0x0037 ) )
    {
    str.str("");
    str << m_Direction[0][0] << "\\"
      << m_Direction[0][1] << "\\";
    /*
     * This is where the 3rd component of the direction is being lost
     * ITK mechanism does not support 2D image, placed in 3D world...
     */
    if( m_Direction.size() == 3 )
      {
      str << m_Direction[0][2] << "\\";
      }
    else
      {
      str << 0. << "\\";
      }
    str << m_Direction[1][0] << "\\"
      << m_Direction[1][1] << "\\";
    if( m_Direction.size() == 3 )
      {
      str << m_Direction[1][2];
      }
    else
      {
      str << 0.;
      }
    header->InsertValEntry(str.str(),0x0020,0x0037); // Image Orientation (Patient)
    }

=======================================================================

Moreover, you could have a look at itk::ImageFileReader::GenerateOutputInformation() and itk::ImageFileWriter::Write() to see that the former inverts the matrix formed by itk::ImageIOBase::m_Direction when setting the image direction of the output image (which is totally right) and the latter does exactly the opposite (which also is totally right). Because of that, itk::GDCMImageIO::Write() may not consider the direction cosines to be stored as columns of the matrix itk::ImageIOBase::m_Direction! This would be true for itk::ImageBase::m_Direction but not for itk::ImageIOBase::m_Direction!!!

Hope you agree with me (and in fact hoping that I'm right :-) ).

Thankx for your efforts!
Andreas

(0015849)
Andreas Schuh (reporter)
2009-03-30 09:28

The short hand explanation should be that the member m_Direction of itk::GDCMImageIO (inherited from itk::ImageIOBase) is a vector that points to vectors, where the first vector contains the direction cosines of the image rows, the second the direction cosines of the image columns and the third contains the direction cosines of the image slices. So, m_Direction[0][0], m_Direction[0][1], m_Direction[0][2] are the direction cosines of the image rows. itk::GDCMImageIO::m_Direction stores the direction cosines as rows.

Now having a look at itk::GDCMImageIO::Write() it turns out, that m_Direction[0][0], m_Direction[1][0], m_Direction[2][0] is taken as direction cosines of the image rows instead ... this mistake might just been made as the itk::ImageBase::m_Direction matrix stores the direction cosines as columns. But the itk::ImageIOBase doesn't!

--
regards
Andreas

(0017624)
Mathieu Malaterre (developer)
2009-09-19 06:49

The attached example will not work in gdcm2 since Secondary Capture Image Storage will be written. Instead one need to say that Modality is -say- MR.

For instance:

    // Be sure that no meta information is present
    ImageType::Pointer image = reader1->GetOutput();
    itk::MetaDataDictionary dict;

  std::string tagkey = "0008|0060"; // Modality
  std::string value = "MR";
  itk::EncapsulateMetaData<std::string>(dict, tagkey, value );

    image->SetMetaDataDictionary( dict );

    // Write image as DICOM file
    WriterType::Pointer writer = WriterType::New();

    writer->SetFileName( argv[2] );
    writer->SetInput( image );

(0017625)
Mathieu Malaterre (developer)
2009-09-19 06:50

$ cvs ci -m"BUG: Fix 0007748: itk::GDCMImageIO::Write() Image Orientation (Patient)" itkGDCMImageIO.cxx ~/Projects/Insight/Code/IO
Committer: Mathieu Malaterre <mathieu.malaterre@gmail.com>
Running style check
Running style check
Processing itkGDCMImageIO.cxx

/cvsroot/Insight/Insight/Code/IO/itkGDCMImageIO.cxx,v <-- itkGDCMImageIO.cxx
new revision: 1.162; previous revision: 1.161

(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.

(0018657)
Luis Ibanez (manager)
2009-12-02 15:10

The fix has now been back ported to the ITK 3.16 branch,
details are in Bug 9993.

Issue History
Date Modified	Username	Field	Change
2008-09-27 19:02	Andreas Schuh	New Issue
2008-10-06 05:31	Mathieu Malaterre	Status	new => assigned
2008-10-06 05:31	Mathieu Malaterre	Assigned To	=> Mathieu Malaterre
2008-11-30 08:56	Mathieu Malaterre	Note Added: 0014233
2008-12-01 18:12	Andreas Schuh	File Added: GDCMImageReadWriteDemo.cxx
2008-12-01 18:29	Andreas Schuh	Note Added: 0014264
2008-12-01 18:38	Andreas Schuh	Note Added: 0014265
2009-02-10 18:58	Andreas Schuh	Note Added: 0014889
2009-03-30 09:28	Andreas Schuh	Note Added: 0015849
2009-09-19 06:49	Mathieu Malaterre	Note Added: 0017624
2009-09-19 06:50	Mathieu Malaterre	Note Added: 0017625
2009-09-19 06:51	Mathieu Malaterre	Note Added: 0017626
2009-09-19 06:51	Mathieu Malaterre	Status	assigned => closed
2009-09-19 06:51	Mathieu Malaterre	Resolution	open => fixed
2009-12-02 15:10	Luis Ibanez	Relationship added	parent of 0009993
2009-12-02 15:10	Luis Ibanez	Note Added: 0018657
2009-12-03 06:27	Mathieu Malaterre	Relationship added	related to 0009622

Notes
(0014233) Mathieu Malaterre (developer) 2008-11-30 08:56	Please attach demo, thanks.

(0014264) Andreas Schuh (reporter) 2008-12-01 18:29	The Demo is quite simple. Try it with a volumetric DICOM file with the DICOM entry (0020,0037) Image Orientation (Patient) be 0/0/-1/1/0/0 for instance. GDCMImageIO::InternalReadImageInformation() reads in the patient orientation correctly and sets it with ImageIOBase::SetDirection() for all three axes. The ImageFileReader then inverts this direction matrix by just setting the row vectors of ImageIOBase::m_Direction as columns of the image direction (which is the inverse patient orientation, the mapping from image space to physical space). Accordingly, ImageFileWriter takes the direction matrix of the image (ImageBase::GetDirection()) and sets the columns of this matrix as rows of the m_Direction member of the ImageIOBase subclass (in our case GDCMImageIO) via ImageIOBase::SetDirection(). Till that point everything is fine. But then in GDCMImageIO::Write() the following code writes the wrong elements of m_Direction as values of the (0020,0037) header entry: if( !header->GetValEntry(0x0020,0x0037 ) ) { str.str(""); str << m_Direction[0][0] << "\\" << m_Direction[1][0] << "\\"; /* * This is where the 3rd component of the direction is being lost * ITK mechanism does not support 2D image, placed in 3D world... */ if( m_Direction.size() == 3 ) { str << m_Direction[2][0] << "\\"; } else { str << 0. << "\\"; } str << m_Direction[0][1] << "\\" << m_Direction[1][1] << "\\"; if( m_Direction.size() == 3 ) { str << m_Direction[2][1]; } else { str << 0.; } header->InsertValEntry(str.str(),0x0020,0x0037); } Actually, (m_Direction[0][0], m_Direction[0][1], m_Direction[0][2]) is the direction of the first axis and (m_Direction[1][0], m_Direction[1][1], m_Direction[1][2]) the one of the second axis. Exactly how ImageFileWriter sets the directions via ImageIOBase::SetDirection().

(0014265) Andreas Schuh (reporter) 2008-12-01 18:38	To resolve some confusion: Patient Orientation: 0/0/-1/1/0/0 GDCMImageIO::m_Direction: 0 0 -1 1 0 0 0 -1 0 Image::m_Direction: 0 1 0 0 0 -1 -1 0 0 The conversion/inversion from one to another is done by ImageFileReader, ImageFileWriter.

(0014889) Andreas Schuh (reporter) 2009-02-10 18:58	Hi Mathieu, maybe you still don't believe me :-) ... Let me explain it another way ... itk::GDCMImageIO::InternalReadImageInformation() sets the m_Direction member of the itk::ImageIOBase class doing the following: m_Direction.resize(3); this->SetDirection(0, rowDirection); this->SetDirection(1, columnDirection); this->SetDirection(2, sliceDirection); So far so good. Having a look at itk::ImageIOBase::SetDirection() reveals the following: m_Direction[0] stores the rowDirection m_Direction[1] stores the columnDirection m_Direction[2] stores the sliceDirection So, rowDirection = [m_Direction[0][0], m_Direction[0][1], m_Direction[0][2]] columnDirection = [m_Direction[1][0], m_Direction[1][1], m_Direction[1][2]] sliceDirection = [m_Direction[2][0], m_Direction[2][1], m_Direction[2][2]] Any objections? Don't think so. So proceeding to the BUGGY code in itk::GDCMImageIO::Write(). ///////////////////////////////// // ATTENTION: THIS IS WRONG!!! // ///////////////////////////////// // Handle Direction = Image Orientation Patient // Same comment as above, if user tell us what the Orientation is, we should not try // to set if from the Image as we might have lost some information if( !header->GetValEntry(0x0020,0x0037 ) ) { str.str(""); str << m_Direction[0][0] << "\\" << m_Direction[1][0] << "\\"; /* * This is where the 3rd component of the direction is being lost * ITK mechanism does not support 2D image, placed in 3D world... / if( m_Direction.size() == 3 ) { str << m_Direction[2][0] << "\\"; } else { str << 0. << "\\"; } str << m_Direction[0][1] << "\\" << m_Direction[1][1] << "\\"; if( m_Direction.size() == 3 ) { str << m_Direction[2][1]; } else { str << 0.; } header->InsertValEntry(str.str(),0x0020,0x0037); // Image Orientation (Patient) } As you see, [m_Direction[0][0], m_Direction[1][0], m_Direction[2][0]] is written as row direction and [m_Direction[0][1], m_Direction[1][1], m_Direction[2][1]] as column direction ... but this is NOT right as we just saw by looking at itk::GDCMImageIO::InternalReadImageInformation(). Instead the code in itk::GDCMImageIO::Write() should look like this: /////////////////////// // HERE'S THE BUGFIX // /////////////////////// // Handle Direction = Image Orientation Patient // Same comment as above, if user tell us what the Orientation is, we should not try // to set if from the Image as we might have lost some information if( !header->GetValEntry(0x0020,0x0037 ) ) { str.str(""); str << m_Direction[0][0] << "\\" << m_Direction[0][1] << "\\"; / * This is where the 3rd component of the direction is being lost * ITK mechanism does not support 2D image, placed in 3D world... */ if( m_Direction.size() == 3 ) { str << m_Direction[0][2] << "\\"; } else { str << 0. << "\\"; } str << m_Direction[1][0] << "\\" << m_Direction[1][1] << "\\"; if( m_Direction.size() == 3 ) { str << m_Direction[1][2]; } else { str << 0.; } header->InsertValEntry(str.str(),0x0020,0x0037); // Image Orientation (Patient) } ======================================================================= Moreover, you could have a look at itk::ImageFileReader::GenerateOutputInformation() and itk::ImageFileWriter::Write() to see that the former inverts the matrix formed by itk::ImageIOBase::m_Direction when setting the image direction of the output image (which is totally right) and the latter does exactly the opposite (which also is totally right). Because of that, itk::GDCMImageIO::Write() may not consider the direction cosines to be stored as columns of the matrix itk::ImageIOBase::m_Direction! This would be true for itk::ImageBase::m_Direction but not for itk::ImageIOBase::m_Direction!!! Hope you agree with me (and in fact hoping that I'm right :-) ). Thankx for your efforts! Andreas

(0015849) Andreas Schuh (reporter) 2009-03-30 09:28	The short hand explanation should be that the member m_Direction of itk::GDCMImageIO (inherited from itk::ImageIOBase) is a vector that points to vectors, where the first vector contains the direction cosines of the image rows, the second the direction cosines of the image columns and the third contains the direction cosines of the image slices. So, m_Direction[0][0], m_Direction[0][1], m_Direction[0][2] are the direction cosines of the image rows. itk::GDCMImageIO::m_Direction stores the direction cosines as rows. Now having a look at itk::GDCMImageIO::Write() it turns out, that m_Direction[0][0], m_Direction[1][0], m_Direction[2][0] is taken as direction cosines of the image rows instead ... this mistake might just been made as the itk::ImageBase::m_Direction matrix stores the direction cosines as columns. But the itk::ImageIOBase doesn't! -- regards Andreas

(0017624) Mathieu Malaterre (developer) 2009-09-19 06:49	The attached example will not work in gdcm2 since Secondary Capture Image Storage will be written. Instead one need to say that Modality is -say- MR. For instance: // Be sure that no meta information is present ImageType::Pointer image = reader1->GetOutput(); itk::MetaDataDictionary dict; std::string tagkey = "0008\|0060"; // Modality std::string value = "MR"; itk::EncapsulateMetaData<std::string>(dict, tagkey, value ); image->SetMetaDataDictionary( dict ); // Write image as DICOM file WriterType::Pointer writer = WriterType::New(); writer->SetFileName( argv[2] ); writer->SetInput( image );

(0017625) Mathieu Malaterre (developer) 2009-09-19 06:50	$ cvs ci -m"BUG: Fix 0007748: itk::GDCMImageIO::Write() Image Orientation (Patient)" itkGDCMImageIO.cxx ~/Projects/Insight/Code/IO Committer: Mathieu Malaterre <mathieu.malaterre@gmail.com> Running style check Running style check Processing itkGDCMImageIO.cxx /cvsroot/Insight/Insight/Code/IO/itkGDCMImageIO.cxx,v <-- itkGDCMImageIO.cxx new revision: 1.162; previous revision: 1.161

(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.

(0018657) Luis Ibanez (manager) 2009-12-02 15:10	The fix has now been back ported to the ITK 3.16 branch, details are in Bug 9993.