Correct coordinate transformations for DICOM data set.

classic Classic list List threaded Threaded
8 messages Options
Reply | Threaded
Open this post in threaded view
|

Correct coordinate transformations for DICOM data set.

Vipul Pai Raikar
Coordinate transformations for DICOM data set.

Hello everybody,

Before I begin to explain my roadblock\problem, let me briefly put forward what I am trying to achieve. I am using VTK to handle the visualization in my work with Image-Guidance (IGS). At the current time I am trying to register a tracked tool (already done) to a data-set and render it correctly in 3D space on screen.

To begin with I have extensively gone through the posts here on understanding how to visualize the 3D space in DICOM coordinate frame (Patient\Scanner\Reference).  To begin with, I am generating a vtkMatrix4X4 with direction cosines and the origin from the dicom header of the first image in the dataset (LowerLeft is ON). I have then applied this as a user matrix to my vtkVolume, which according to my understanding from all the wonderful explanations here (Especially Dr Gobbi), will apply the DICOM coordinate frame to the rendering (for the lack of better way to word this).

My questions are two fold. First: How to I now render a point, which is registered to this coordinate frame, correctly (in the form of a vtkActor)?
Second: How to I move this object\actor in real-time?

I have tried and failed to correctly do so by applying the same user transform to the actor. I am still quite new to this. I would very much appreciate any help I can get to move forward in the right direction.

best regards,

V
Reply | Threaded
Open this post in threaded view
|

Re: Correct coordinate transformations for DICOM data set.

David Gobbi
Hi V,

I found the vtkDICOMImageReader to be so lacking in functionality
that I wrote my own: http://dgobbi.github.io/vtk-dicom/

And, of course, many other people use ITK to read their DICOMs
(though you have to be very careful or you will lose coord system
fidelity when importing the image from ITK to VTK).  Finally, there
is also the vtkGDCMImageReader.


But if you are using the vtkDICOMImageReader, here are a few tips:

1) The reader ignores FileLowerLeft, so turning it on does nothing

2) The reader always applies a vertical flip to the images (nearly all
VTK image readers flip images when they read them).

3) The reader uses the "Image Position (Patient)" to sort the slices,
but it reverses the order.  It has to do this, because it has already
done a "y" flip, and two flips in total are necessary to avoid having
the volume become a mirror image of itself.

4) When you call reader->GetImagePositionPatient(), it returns the
information from the last slice that it read.  Since it reads the
slices in reverse order (see #3), this actually means that it returns
the information from the first slice.


So, putting everything together, here is the code that I used to read
images with the vtkDICOMImageReader (before I wrote my own reader):

  vtkSmartPointer<vtkDICOMImageReader> reader =
    vtkSmartPointer<vtkDICOMImageReader>::New();
  reader->SetDirectoryName(directoryName);

  // flip the image in Y and Z directions
  vtkSmartPointer<vtkImageReslice> flip =
    vtkSmartPointer<vtkImageReslice>::New();
  flip->SetInputConnection(reader->GetOutputPort());
  flip->SetResliceAxesDirectionCosines(
    1,0,0, 0,-1,0, 0,0,-1);

  flip->Update();
  vtkSmartPointer<vtkImageData> image = flip->GetOutput();

  // generate the matrix
  float *position = reader->GetImagePositionPatient();
  float *orientation = reader->GetImageOrientationPatient();
  float *xdir = &orientation[0];
  float *ydir = &orientation[3];
  float zdir[3];
  vtkMath::Cross(xdir, ydir, zdir);

  vtkSmartPointer<vtkMatrix4x4> matrix =
    vtkSmartPointer<vtkMatrix4x4>::New();
  for (int i = 0; i < 3; i++)
    {
    matrix->Element[i][0] = xdir[i];
    matrix->Element[i][1] = ydir[i];
    matrix->Element[i][2] = zdir[i];
    matrix->Element[i][3] = position[i];
    }


Then you can attach the matrix to the vtkVolume with
SetUserMatrix(matrix).  If you have points that are already in
the DICOM patient coordinate system, then you don't need to
apply any transformations to those points.  I'm not going to give
detailed info about how to display points in VTK (you can probably
find wiki examples), but my advice is to use vtkGlyph3D with
vtkSphereSource.  Just search for "VTK Glyph3D" and you'll
find a suitable example.

 - David


On Wed, Aug 20, 2014 at 9:55 AM, vipulraikar <[hidden email]> wrote:

> Coordinate transformations for DICOM data set.
>
> Hello everybody,
>
> Before I begin to explain my roadblock\problem, let me briefly put forward
> what I am trying to achieve. I am using VTK to handle the visualization in
> my work with Image-Guidance (IGS). At the current time I am trying to
> register a tracked tool (already done) to a data-set and render it correctly
> in 3D space on screen.
>
> To begin with I have extensively gone through the posts here on
> understanding how to visualize the 3D space in DICOM coordinate frame
> (Patient\Scanner\Reference).  To begin with, I am generating a vtkMatrix4X4
> with direction cosines and the origin from the dicom header of the first
> image in the dataset (LowerLeft is ON). I have then applied this as a user
> matrix to my vtkVolume, which according to my understanding from all the
> wonderful explanations here (Especially Dr Gobbi), will apply the DICOM
> coordinate frame to the rendering (for the lack of better way to word this).
>
> My questions are two fold. First: How to I now render a point, which is
> registered to this coordinate frame, correctly (in the form of a vtkActor)?
> Second: How to I move this object\actor in real-time?
>
> I have tried and failed to correctly do so by applying the same user
> transform to the actor. I am still quite new to this. I would very much
> appreciate any help I can get to move forward in the right direction.
>
> best regards,
>
> V
_______________________________________________
Powered by www.kitware.com

Visit other Kitware open-source projects at http://www.kitware.com/opensource/opensource.html

Please keep messages on-topic and check the VTK FAQ at: http://www.vtk.org/Wiki/VTK_FAQ

Follow this link to subscribe/unsubscribe:
http://public.kitware.com/mailman/listinfo/vtkusers
Reply | Threaded
Open this post in threaded view
|

Re: Correct coordinate transformations for DICOM data set.

David Gobbi
About the code that I included in my previous email, I was
sure that I would forget something, and I did!

5) The "Origin" of the vtkImageData that you get from the
DICOM reader _must_ be (0,0,0) or else the matrix will not
work.  The matrix already has the position information from
the DICOM header.  The VTK "Origin" is not the same as the
DICOM "Position", so you cannot get "Image Position (Patient)"
from the DICOM header and then use that as your Origin!

Likewise, if you ever import oriented ITK images into VTK,
the ITK Origin and the VTK Origin are not the same thing!
The ITK Origin is defined with respect to the oriented image
(I.e. with respect to a rotated coordinate system), while the
VTK Origin is defined with respect to an un-oriented image
(i.e. with respect to a coordinate system where X, Y, and Z
are aligned with the column, row, and slice directions).

I hope this clarifies things, rather than just confusing them
further!

 - David



On Wed, Aug 20, 2014 at 10:58 AM, David Gobbi <[hidden email]> wrote:

> Hi V,
>
> I found the vtkDICOMImageReader to be so lacking in functionality
> that I wrote my own: http://dgobbi.github.io/vtk-dicom/
>
> And, of course, many other people use ITK to read their DICOMs
> (though you have to be very careful or you will lose coord system
> fidelity when importing the image from ITK to VTK).  Finally, there
> is also the vtkGDCMImageReader.
>
>
> But if you are using the vtkDICOMImageReader, here are a few tips:
>
> 1) The reader ignores FileLowerLeft, so turning it on does nothing
>
> 2) The reader always applies a vertical flip to the images (nearly all
> VTK image readers flip images when they read them).
>
> 3) The reader uses the "Image Position (Patient)" to sort the slices,
> but it reverses the order.  It has to do this, because it has already
> done a "y" flip, and two flips in total are necessary to avoid having
> the volume become a mirror image of itself.
>
> 4) When you call reader->GetImagePositionPatient(), it returns the
> information from the last slice that it read.  Since it reads the
> slices in reverse order (see #3), this actually means that it returns
> the information from the first slice.
>
>
> So, putting everything together, here is the code that I used to read
> images with the vtkDICOMImageReader (before I wrote my own reader):
>
>   vtkSmartPointer<vtkDICOMImageReader> reader =
>     vtkSmartPointer<vtkDICOMImageReader>::New();
>   reader->SetDirectoryName(directoryName);
>
>   // flip the image in Y and Z directions
>   vtkSmartPointer<vtkImageReslice> flip =
>     vtkSmartPointer<vtkImageReslice>::New();
>   flip->SetInputConnection(reader->GetOutputPort());
>   flip->SetResliceAxesDirectionCosines(
>     1,0,0, 0,-1,0, 0,0,-1);
>
>   flip->Update();
>   vtkSmartPointer<vtkImageData> image = flip->GetOutput();
>
>   // generate the matrix
>   float *position = reader->GetImagePositionPatient();
>   float *orientation = reader->GetImageOrientationPatient();
>   float *xdir = &orientation[0];
>   float *ydir = &orientation[3];
>   float zdir[3];
>   vtkMath::Cross(xdir, ydir, zdir);
>
>   vtkSmartPointer<vtkMatrix4x4> matrix =
>     vtkSmartPointer<vtkMatrix4x4>::New();
>   for (int i = 0; i < 3; i++)
>     {
>     matrix->Element[i][0] = xdir[i];
>     matrix->Element[i][1] = ydir[i];
>     matrix->Element[i][2] = zdir[i];
>     matrix->Element[i][3] = position[i];
>     }
>
>
> Then you can attach the matrix to the vtkVolume with
> SetUserMatrix(matrix).  If you have points that are already in
> the DICOM patient coordinate system, then you don't need to
> apply any transformations to those points.  I'm not going to give
> detailed info about how to display points in VTK (you can probably
> find wiki examples), but my advice is to use vtkGlyph3D with
> vtkSphereSource.  Just search for "VTK Glyph3D" and you'll
> find a suitable example.
>
>  - David
>
>
> On Wed, Aug 20, 2014 at 9:55 AM, vipulraikar <[hidden email]> wrote:
>> Coordinate transformations for DICOM data set.
>>
>> Hello everybody,
>>
>> Before I begin to explain my roadblock\problem, let me briefly put forward
>> what I am trying to achieve. I am using VTK to handle the visualization in
>> my work with Image-Guidance (IGS). At the current time I am trying to
>> register a tracked tool (already done) to a data-set and render it correctly
>> in 3D space on screen.
>>
>> To begin with I have extensively gone through the posts here on
>> understanding how to visualize the 3D space in DICOM coordinate frame
>> (Patient\Scanner\Reference).  To begin with, I am generating a vtkMatrix4X4
>> with direction cosines and the origin from the dicom header of the first
>> image in the dataset (LowerLeft is ON). I have then applied this as a user
>> matrix to my vtkVolume, which according to my understanding from all the
>> wonderful explanations here (Especially Dr Gobbi), will apply the DICOM
>> coordinate frame to the rendering (for the lack of better way to word this).
>>
>> My questions are two fold. First: How to I now render a point, which is
>> registered to this coordinate frame, correctly (in the form of a vtkActor)?
>> Second: How to I move this object\actor in real-time?
>>
>> I have tried and failed to correctly do so by applying the same user
>> transform to the actor. I am still quite new to this. I would very much
>> appreciate any help I can get to move forward in the right direction.
>>
>> best regards,
>>
>> V
_______________________________________________
Powered by www.kitware.com

Visit other Kitware open-source projects at http://www.kitware.com/opensource/opensource.html

Please keep messages on-topic and check the VTK FAQ at: http://www.vtk.org/Wiki/VTK_FAQ

Follow this link to subscribe/unsubscribe:
http://public.kitware.com/mailman/listinfo/vtkusers
Reply | Threaded
Open this post in threaded view
|

Re: Correct coordinate transformations for DICOM data set.

Vipul Pai Raikar
Thank you for your prompt and detailed response Dr Gobbi. I was just about to ask you about the origin. I shall try and implement these ideas and get back.

Thank you very much!

V
Reply | Threaded
Open this post in threaded view
|

Re: Correct coordinate transformations for DICOM data set.

David Gobbi
Just to be complete, here is a code block that includes a filter
that forces the Origin to be (0,0,0), so that all positioning is done
by the matrix:

  vtkSmartPointer<vtkDICOMImageReader> reader =
    vtkSmartPointer<vtkDICOMImageReader>::New();
  reader->SetDirectoryName(directoryName);

  // flip the image in Y and Z directions
  vtkSmartPointer<vtkImageReslice> flip =
    vtkSmartPointer<vtkImageReslice>::New();
  flip->SetInputConnection(reader->GetOutputPort());
  flip->SetResliceAxesDirectionCosines(
    1,0,0, 0,-1,0, 0,0,-1);

  // force the origin to be (0,0,0) since we set position in the matrix
  vtkSmartPointer<vtkImageChangeInformation> change =
    vtkSmartPointer<vtkImageChangeInformation>::New();
  change->SetInputConnection(flip->GetOutputPort());
  change->SetOutputOrigin(0,0,0);

  change->Update();
  vtkSmartPointer<vtkImageData> image = change->GetOutput();

  // generate the matrix
  float *position = reader->GetImagePositionPatient();
  float *orientation = reader->GetImageOrientationPatient();
  float *xdir = &orientation[0];
  float *ydir = &orientation[3];
  float zdir[3];
  vtkMath::Cross(xdir, ydir, zdir);

  vtkSmartPointer<vtkMatrix4x4> matrix =
    vtkSmartPointer<vtkMatrix4x4>::New();
  for (int i = 0; i < 3; i++)
    {
    matrix->Element[i][0] = xdir[i];
    matrix->Element[i][1] = ydir[i];
    matrix->Element[i][2] = zdir[i];
    matrix->Element[i][3] = position[i];
    }

Definitely this is a lot of work just to read a DICOM image in the
correct coordinate system!  Let me know if this gives you the
results that you need.

 - David

On Wed, Aug 20, 2014 at 11:21 AM, Vipul Pai Raikar <[hidden email]> wrote:

> Thank you for your prompt and detailed response Dr Gobbi. I was just about to
> ask you about the origin. I shall try and implement these ideas and get
> back.
>
> Thank you very much!
>
> V
>
>
>
> --
> View this message in context: http://vtk.1045678.n5.nabble.com/Correct-coordinate-transformations-for-DICOM-data-set-tp5728309p5728313.html
> Sent from the VTK - Users mailing list archive at Nabble.com.
> _______________________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at http://www.kitware.com/opensource/opensource.html
>
> Please keep messages on-topic and check the VTK FAQ at: http://www.vtk.org/Wiki/VTK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://public.kitware.com/mailman/listinfo/vtkusers
_______________________________________________
Powered by www.kitware.com

Visit other Kitware open-source projects at http://www.kitware.com/opensource/opensource.html

Please keep messages on-topic and check the VTK FAQ at: http://www.vtk.org/Wiki/VTK_FAQ

Follow this link to subscribe/unsubscribe:
http://public.kitware.com/mailman/listinfo/vtkusers
Reply | Threaded
Open this post in threaded view
|

Re: Correct coordinate transformations for DICOM data set.

Vipul Pai Raikar
Hello Dr Gobbi,

Well as you explained, I implemented the solution. It works like a charm. I agree with you that its quite a few lines of code to get it to work, but it indeed helped my better grasp the concept.
As well as it worked, there was however an issue with the rendering being flipped in the Y-direction. I approached this problem by manipulating the camera associated with the renderer. I remember from some of your earlier posts that I read here, you explicitly advise setting up the camera by your self. As a result the rendering is now oriented correctly. I was just wondering if this was a normal thing or is the flip the byproduct of something else? ; and if this approach is correct.
I also can now correctly position my tracker data.

best regards,

V
Reply | Threaded
Open this post in threaded view
|

Re: Correct coordinate transformations for DICOM data set.

David Gobbi
Hi Vipul,

For the default position of the VTK camera,

X points right
Y points up
Z points out of the screen

So, yes, if one simply tries to slap a DICOM image up on
the screen with VTK, it will be displayed with the wrong
orientation.  After all, VTK is a general-purpose graphics
library, so its default settings are not tuned specifically
for DICOM.

Setting up the camera properly, so that you get the desired
view of your data, is a must.

 - David

On Thu, Aug 21, 2014 at 11:14 AM, Vipul Pai Raikar <[hidden email]> wrote:

> Hello Dr Gobbi,
>
> Well as you explained, I implemented the solution. It works like a charm. I
> agree with you that its quite a few lines of code to get it to work, but it
> indeed helped my better grasp the concept.
> As well as it worked, there was however an issue with the rendering being
> flipped in the Y-direction. I approached this problem by manipulating the
> camera associated with the renderer. I remember from some of your earlier
> posts that I read here, you explicitly advise setting up the camera by your
> self. As a result the rendering is now oriented correctly. I was just
> wondering if this was a normal thing or is the flip the byproduct of
> something else? ; and if this approach is correct.
> I also can now correctly position my tracker data.
>
> best regards,
>
> V
>
>
>
> --
> View this message in context: http://vtk.1045678.n5.nabble.com/Correct-coordinate-transformations-for-DICOM-data-set-tp5728309p5728324.html
> Sent from the VTK - Users mailing list archive at Nabble.com.
> _______________________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at http://www.kitware.com/opensource/opensource.html
>
> Please keep messages on-topic and check the VTK FAQ at: http://www.vtk.org/Wiki/VTK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://public.kitware.com/mailman/listinfo/vtkusers
_______________________________________________
Powered by www.kitware.com

Visit other Kitware open-source projects at http://www.kitware.com/opensource/opensource.html

Please keep messages on-topic and check the VTK FAQ at: http://www.vtk.org/Wiki/VTK_FAQ

Follow this link to subscribe/unsubscribe:
http://public.kitware.com/mailman/listinfo/vtkusers
Reply | Threaded
Open this post in threaded view
|

Re: Correct coordinate transformations for DICOM data set.

Vipul Pai Raikar
Thank you very much for all the explanations. Setting up the camera correctly has done the trick.

Best regards,

V