How to get voxel coordinates(as in ITK snap ) from vtkimage viewer from mouse click

Engineering
, ,
1

Hello,
I am trying to use ITK’s region growing algorithm for medical image segmentation. So to put the seedpoint i need exact voxel coordinates of the mouse click position in my vtkimageviewer. So how can I get those voxel coordinates(as in ITK snap tool) from my vtkimageviewer. Currently I am able to get VTK World Coordinates using pick point picker . Now how to convert this vtk world coordinates to voxel coordinates as in ITK snap tool …
Please help me with this information.

2

Have you tried my suggestion at stack overflow? If so, why doesn’t that work?

There was an example which I wrote:
https://itk.org/Wiki/ITK/Examples/IO/itkVtkImageConvertDICOM
It has been since deleted as unmaintained. Here is the code:

#include <itkRandomImageSource.h>
#include "itkImageToVTKImageFilter.h"
 
#include "vtkVersion.h"
#include <vtkSmartPointer.h>
#include <vtkGPUVolumeRayCastMapper.h>
#include <vtkColorTransferFunction.h>
#include <vtkPiecewiseFunction.h>
#include <vtkRenderer.h>
#include <vtkRenderWindow.h>
#include <vtkRenderWindowInteractor.h>
#include <vtkVolumeProperty.h>
#include <vtkMatrix4x4.h>
#include <vtkAxesActor.h>

using VisualizingImageType = itk::Image<unsigned char, 3>;

VisualizingImageType::Pointer createImage()
{
    itk::Size<3> size;
    size[0]=512;
    size[1]=512;
    size[2]=16;
 
    itk::RandomImageSource<VisualizingImageType>::Pointer randomImageSource =
        itk::RandomImageSource<VisualizingImageType>::New();
    randomImageSource->SetNumberOfThreads(1); // to produce non-random results
    randomImageSource->SetSize(size);
    randomImageSource->Update();
    
    VisualizingImageType::Pointer image=randomImageSource->GetOutput();
    double t[3]={0, -200, 150};
    image->SetOrigin(t);
    t[0]=0.68;
    t[1]=0.68;
    t[2]=4.4;
    image->SetSpacing(t);

    VisualizingImageType::DirectionType m;
    m(0,0)=0;
    m(0,1)=1;
    m(0,2)=0;
    m(1,0)=0.707;
    m(1,1)=0;
    m(1,2)=-0.707;
    m(2,0)=-0.707;
    m(2,1)=0;
    m(2,2)=-0.707;
    m=m.GetTranspose();
    image->SetDirection(m);
    return image;
}

int main(int argc, char *argv[])
{
  VisualizingImageType::Pointer image=createImage();
 
  vtkSmartPointer<vtkRenderWindow> renWin = vtkSmartPointer<vtkRenderWindow>::New();
  vtkSmartPointer<vtkRenderer> ren1 = vtkSmartPointer<vtkRenderer>::New();
  ren1->SetBackground(0.5f,0.5f,1.0f);
 
  renWin->AddRenderer(ren1);
  renWin->SetSize(1280,1024);
  vtkSmartPointer<vtkRenderWindowInteractor> iren =
      vtkSmartPointer<vtkRenderWindowInteractor>::New();
  iren->SetRenderWindow(renWin);
  renWin->Render(); // make sure we have an OpenGL context.
 
  using itkVtkConverter = itk::ImageToVTKImageFilter<VisualizingImageType>;
  itkVtkConverter::Pointer conv=itkVtkConverter::New();
  conv->SetInput(image);
 
  vtkSmartPointer<vtkGPUVolumeRayCastMapper> volumeMapper =
      vtkSmartPointer<vtkGPUVolumeRayCastMapper>::New();
#if VTK_MAJOR_VERSION <= 5
  volumeMapper->SetInput(conv->GetOutput());
#else
  conv->Update();
  volumeMapper->SetInputData(conv->GetOutput());
#endif
 
  vtkSmartPointer<vtkVolumeProperty> volumeProperty =
  vtkSmartPointer<vtkVolumeProperty>::New();
 
  vtkSmartPointer<vtkPiecewiseFunction> compositeOpacity =
  vtkSmartPointer<vtkPiecewiseFunction>::New();
  compositeOpacity->AddPoint(0.0,0.1);
  compositeOpacity->AddPoint(80.0,0.2);
  compositeOpacity->AddPoint(255.0,0.1);
  volumeProperty->SetScalarOpacity(compositeOpacity);
 
  vtkSmartPointer<vtkColorTransferFunction> color =
  vtkSmartPointer<vtkColorTransferFunction>::New();
  color->AddRGBPoint(0.0  ,0.0,0.0,1.0);
  color->AddRGBPoint(40.0  ,1.0,0.0,0.0);
  color->AddRGBPoint(255.0,1.0,1.0,1.0);
  volumeProperty->SetColor(color);
 
  vtkSmartPointer<vtkVolume> volume =
  vtkSmartPointer<vtkVolume>::New();
  volume->SetMapper(volumeMapper);
  volume->SetProperty(volumeProperty);
 
  // Here we take care of position and orientation
  // so that volume is in DICOM patient physical space
  VisualizingImageType::DirectionType d=image->GetDirection();
  vtkMatrix4x4 *mat=vtkMatrix4x4::New(); //start with identity matrix
  for (int i=0; i<3; i++)
      for (int k=0; k<3; k++)
          mat->SetElement(i,k, d(i,k));
 
  // Counteract the built-in translation by origin
  VisualizingImageType::PointType origin=image->GetOrigin();
  volume->SetPosition(-origin[0], -origin[1], -origin[2]);
 
  // Add translation to the user matrix
  for (int i=0; i<3; i++)
    {
    mat->SetElement(i,3, origin[i]);
    }
  volume->SetUserMatrix(mat);
 
  // Add coordinate system axes, so we have a reference for position and orientation
  vtkSmartPointer<vtkAxesActor> axes = vtkSmartPointer<vtkAxesActor>::New();
  axes->SetTotalLength(250,250,250);
  axes->SetShaftTypeToCylinder();
  axes->SetCylinderRadius(0.01);
  ren1->AddActor(axes);
 
  ren1->AddVolume( volume );
  ren1->ResetCamera();
 
  renWin->Render();
  renWin->Render();
  renWin->Render();
 
  iren->Start();
 
  return EXIT_SUCCESS;
}

You might need something akin to volume->SetUserMatrix(mat); for ijkIndex = itkImage->TransformPhysicalPointToIndex(worldCoordinates); to work directly. Or an inverse version of that, going from VTK’s world coordinates to ITK’s DICOM Image Position Patient coordinate system before you can use those coordinates in TransformPhysicalPointToIndex.

3

Thank you so much for your response… I will be trying this soon and will give you the update.

4

Hi ,
I am really sorry for the late reply. I tried what you suggested but still it is not giving me the exact coordinates. I don’t know where I am going wrong. What I was doing was,
I am having multiplanar reconstruction window(MPR) and I am using vtkResliceCursorWidget for getting the coordinates . I have tried two methods ,
1)

this->RCW[i]->GetResliceCursorRepresentation()->GetPlaneSource()->GetOrigin(Origin);

this->RCW[i]->GetResliceCursorRepresentation()->GetPlaneSource()->GetPoint1(Point1);

this->RCW[i]->GetResliceCursorRepresentation()->GetPlaneSource()->GetPoint2(Point2);

And then passing these points to TransformPhysicalPointToIndex() .

int evenposition[3];
this->Interactor->GetEventPosition(evenposition);
vtkSmartPointer vtkCellPicker pickerr = vtkSmartPointer vtkCellPicker ::New();
this->Interactor->SetPicker(pickerr);
pickerr->Pick(evenposition[0], evenposition[1], 0, this->Interactor->GetRenderWindow()->GetRenderers()->GetFirstRenderer());
double datacoordinates[3];
pickerr->GetMapperPosition(datacoordinates);

Then passing these points to TransformPhysicalPointToIndex() .

Both of the above methods are not giving me any transformations for getting voxel coordinates as in ITK SNAP.
Could you please tell me where I am going wrong and help me to get the voxel coordinates.

5

I don’t have a ready solution for that. Maybe try looking at source code of ITK-SNAP? This seems to be a part of it. Or look at source of Slicer or MITK? Or use one of those for UI, so you don’t have to make your own GUI?

6

Hello ,
Could you please tell is there any other way to solve this issue other than using other GUI. Please see the below code which is my current status.

                    int evenpos[3];
		this->Interactor->GetEventPosition(evenpos);
		vtkSmartPointer<vtkCellPicker> pickerr = vtkSmartPointer<vtkCellPicker>::New();
		this->Interactor->SetPicker(pickerr);
		pickerr->Pick(evenpos[0], evenpos[1], 0, this->GetDefaultRenderer());
		double picked[3];
		pickerr->GetPickPosition(picked);

Now this picked[3] I need to convert to voxel coordinates. Please help me with this

7

Have you tried vtkImageData::TransformPhysicalPointToContinuousIndex()?

If that is not enough, you might need to transform your point picked using the same transforms you apply to your image in VTK scene. Or maybe the inverse of that. Apply transforms before passing the point to TransformPhysicalPointToContinuousIndex.