[vtkusers] vtkDataSet::FindPoint, vtkExtractVOI and world coordinates?

[David Gobbi]

Hi Elvis,

Let's break this down into two parts:

1) Going from VTK world coords to VTK data coords is achieved via the prop
matrix (where a "prop" is a vtkActor, a vtkVolume, a vtkImageSlice, etc):

For data-to-world, apply the matrix to the point.
For world-to-data, apply the inverse of the matrix.

2) Going from VTK data coords (x,y,z) to point indices (i,j,k) is done with
these equations:

x = origin[0] + i * spacing[0]
y = origin[1] + j * spacing[1]
z = origin[2] + k * spacing[2]

i = round((x - origin[0])/spacing[0])
j = round((y - origin[1])/spacing[1])
k = round((z - origin[2])/spacing[2])

Personally, I avoid using any VTK APIs to do these computations, for two
reasons:

a) I prefer mathematical abstractions over programatic abstractions.  If I
have the equations, then I can combine them and simplify them and often
achieve nice performance gains as a result.

b) When I look at the available APIs, sometimes I find they were written
for a subtly different purpose than my own.  For example, consider
vtkImageData::ComputeStructuredPoints:

int ComputeStructuredCoordinates(const double x[3],  int ijk[3],  double
pcoords[3]);

It looks like it might do what you want, right?  But in fact it was written
to work with VTK's abstract concept of a "cell" where, for vtkImageData, a
cell is a box with a sample point (i.e. what we typically think of as a
voxel) at each corner.  So the above method gives us indices to tell us
which cell we are in, but we'd have to check the pcoords to find out what
point we are closest to.  At that point, we might as well just ignore this
API and use the equations given in (2).


As an aside, whenever you're dealing with VTK data, you should always keep
in mind that the APIs were written according to the way that VTK's defines
its "data set" abstraction, i.e. points and cells and attributes thereof.
It's necessary to have a firm grasp of this abstraction in order to use VTK
to its utmost.  However, it's equally necessary to understand that not all
problems fit nicely into this abstraction, and sometimes VTK's API's aren't
going to nicely match the problem that one is trying to solve.  In those
cases, it can be useful to (a) read the VTK source code to understand
exactly what VTK is doing, and (b) figure out exactly how you might have to
do things differently.

 - David



On Wed, Sep 7, 2016 at 9:55 AM, Elvis Stansvik <elvis.stansvik at orexplore.com
> wrote:

> Den 7 sep. 2016 3:32 em skrev "Cory Quammen" <cory.quammen at kitware.com>:
> >
> > Elvis,
> >
> > Have a look at vtkImageData::FindPoint(double x, double y, double z).
> >
> > It takes into account the image spacing and origin, but will not
> > correct for any modifications to the position of the vtkVolume you
> > seem to be using to represent it in the scene.
>
> Thanks Cory,
>
> I did initially find FindPoint, but I was unsure how I can translate the
> vtkIdType it returns (some sort of index?) into i,j,k voxel coordinate. Do
> you know?
>
> Since then I've also found  vtkImageData::ComputeStructuredCoordinates,
> which looks similar, but does seem to produce a ijk voxel coordinate.
> Perhaps that is the function I'm looking for? It's doubly interesting since
> it has a static variant where I can specify the origin, spacing etc
> manually, which I might have a need for actually, if I want to do the same
> computation but don't have access to an actual vtkImageData instance, only
> its parameters so to speak.
>
> Elvis
>
> >
> > HTH,
> > Cory
> >
> >
> >
> > On Wed, Sep 7, 2016 at 3:10 AM, Elvis Stansvik
> > <elvis.stansvik at orexplore.com> wrote:
> > > Den 4 sep. 2016 2:25 em skrev "Elvis Stansvik"
> > > <elvis.stansvik at orexplore.com>:
> > >>
> > >> Hi all,
> > >>
> > >> I have multiple volumes in a renderer, positioned above/below each
> other
> > >> (like a "string" of volumes) along the Z axis, none of them are
> > >> intersecting. As the user clicks+drags the mouse in the viewport, I
> need to
> > >> capture the mouse position at the start/end positions of the drag
> operation
> > >> and convert them to world coordinates. I think I have these parts
> figured
> > >> out.
> > >>
> > >> Then, I need to search through the list of my volumes (sorted by
> their end
> > >> positions in world coordinates), and construct a series of
> vtkExtractVOI
> > >> that will extract the selected segments. For volumes that are
> completely
> > >> inside the selected area, the vtkExtractVOI will extract the full
> volumes,
> > >> but for the start/end volume (which may be the same volume), it will
> only
> > >> extract a portion of it.
> > >>
> > >> I looked at the vtkDataSet::FindPoint, which I think will give me a
> linear
> > >> point ID (not familiar with that from before) for a point given in
> world
> > >> coordinates. But how can I convert that to a structured point
> position that
> > >> I can use when setting up the extent to use for my vtkExtractVOI?
> > >>
> > >> Is the static vtkStructuredData::ComputePointStructuredCoords helper
> the
> > >> function I'm looking for? Or is there some other easier way to
> extract an
> > >> area of a volume when the area is expressed in world coordinates?
> > >> (specifically, as a start/end Z position in world coordinates).
> > >>
> > >> Here's an ASCII drawing that perhaps explains better what I want to
> do,
> > >> the rectangular areas are my volumes:
> > >>
> > >> ^  +-----+
> > >> |  |     |
> > >> |  | P1  |  <--- User clicks here (world coordinate X, Y, Z)..
> > >> |  |     |
> > >> Z  +-----+
> > >>    +-----+
> > >>    |     |
> > >>    |     |
> > >>    +-----+
> > >>    +-----+
> > >>    |     |
> > >>    |     |
> > >>    |  P2 |  <--- ..and drags to here (world coordinate A, B, C).
> > >>    |     |
> > >>    |     |
> > >>    |     |
> > >>    +-----+
> > >>
> > >> From this I want to create three vtkExtractVOI, that when combined
> gives a
> > >> list of vtkImageData that represents the selected area along the Z
> axis (in
> > >> the X and Y direction, everything from the volumes should be
> included).
> > >>
> > >> Thanks in advance,
> > >> Elvis
> > >
> > > I realized my mail above is a little long-winded and unclear. Sorry
> about
> > > that. In essence, what I'm asking is how to map an x,y,x (double)
> point in
> > > world coordinates to a i,j,k (int) grid position (voxel coordinate) of
> a
> > > given volume. I don't want the nearest point, but to return "null" or
> > > equivalent if the world point does not fall inside any voxel.
> > >
> > > I'd also be interested how to do the same but against a plain
> vtkImageData
> > > instead of a vtkVolume. For this I assume I'd have to provide some info
> > > about where in the world the vtkImageData belongs (where/how I intend
> to map
> > > it into the scene).
> > >
> > > Anyone know what the VTK API for this is?
> > >
> > > Thanks in advance,
> > > Elvis
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://public.kitware.com/pipermail/vtkusers/attachments/20160907/9990cbfe/attachment.html>