[Insight-users] How to define and operate on arbitrarily oriented image regions?

[Gavin Baker]

Hello,

Some thoughts on your interesting problems, FWIW...

On Wed, Aug 17, 2005 at 11:56:30AM +0200, Jger Hansegrd wrote:
> I want to define a region in a 3D ITK image. The region is shaped as a
> quadratic prism, but may be oriented arbitrarily in space. That means
> that the sides of the prism are generally not parallel to the image
> axes.
> 
> * Is it possible to define such a region using ITK out of the box? I
> have only been able to find regions that have sides parallel to the
> image axes.

itk::Regions are necessarily tied to the image representation, for obvious
reasons.  If you need this sort of flexibility, you can generate an
itk::Mesh that represents your prism, and use it as a mask.  Then you have
virtually complete flexibility, albeit at the cost of some extra
infrastructure or processing stages.

> On this region, I would like to do the following:
> * Filtering - I am interested in filtering the image only within the
> specified region. Can ITK filters, such as the smoothing filters or
> gradient filters, operate on such regions?

The pipeline infrastructure sets up Filters to operate on an itk::Region,
and AFAIK this is the only means of restricting the processing area in the
pipeline itself.

However, that's not really an issue as you can just filter out the
prism-shaped region of interest.  First, you can calculate your mask and use
the bounding box of the prism to define a Region.  Also, you can resample
the data into your prism and process it entirely.

One idea goes like this... The itk::TriangleMeshToBinaryImageFilter can be
used to render or rasterise a mesh into a volume, which can then be used as
a binary mask (provided you match the image size etc to the other image).
The itk::MaskImageFilter can process an image and produce a masked version
on the output, leaving just the prism-shaped region of the input image.

Something like:

  InputImage                    Mesh (Prism)
      |                          |
      |                TriangleMeshToBinaryImageFilter
      |                          |
      `-------+         +--------'
              |         |
             MaskImageFilter
                    |
                    |
            ResampleImageFilter ------ Transform
                    |
                    |
                    . 
                    .

Where the transform orients the resampled image such that the desired prism
face is "in plane".

> * Sub volume iteration - can any of the image iterators operate only
> on the pixels within this region?

There is quite a variety of iterators, but I don't think there is one that
uses a mesh mask (or SpatialObject mask).

  http://www.itk.org/Doxygen/html/group__ImageIterators.html

However, given a Mesh armed with a method like:

  IsInside(const PointType& p)

it would be quite easy to write an iterator that only operated inside the
mesh.  You could get the bounding box of the mesh to set the requested
region (and speed things up a little), and have some boundary condition to
handle the edges.

Then you'd have a useful general-purpose masking iterator, which we could
then add back into ITK! :)

> * Sub volume data extraction -  I need to extract the sub volume
> included in the region using some type of interpolation, for example
> nearest neigbor or bilinear interpolation. The result should be a new
> image of the same size as the requested region, and with pixel values
> sampled on a regular grid that is oriented along the quadratic prism
> region.

Certainly, the itk::ResampleImageFilter can do this sort of thing; it is
extremely flexible, provided you can express what you want in terms of the
Transforms.  You can resample with the correct orientation, then mask out
just the region you want with the prism.

There are several interpolators available, including Linear, NN, and
BSpline.  The full list is here:

  http://www.itk.org/Doxygen/html/group__ImageInterpolators.html  

> Is this functionality available in ITK, and if so, where should I look
> for documentation?
> If this functionality is not available, I'll have to implement this.
> Any hints on which base classes to use for each type of functionality
> is highly appreciated.

The ImageMaskSpatialObject can also be used in certain contexts, depending
on the algorithm (some discussion appeared on this list in recent weeks).

You could start by creating your own PrismMeshSource derived from
MeshSource.  It could take as parameters the three size dimensions, and an
orientation vector.

For the iterator, there are a few ways you could implement it.  One
possibility is to derive from ConditionalConstIterator, take an input mesh
(which will need its own template parameter), and implement the
IsPixelIncluded() to test if it is inside the mesh.

And yet another way might be to derive from ConditionalConstIterator, take
an ImageSpatialObject as a mask, and then the pixel test just transforms the
pixel into the SO space and returns true if the pixel is non-0.  Then the
image mask can be generated by rasterising the mesh into an image (which is
wrapped by the ImageSpatialObject).  This could be more efficient, as
testing if a point is inside an arbitrary mesh can be expensive, so
rasterising up-front would make this much faster.  Also it means different
types of images could be used as masks.

So something like:


           Image                                   Mesh
              |                                     |
              |                     TriangleMeshToBinaryImageFilter
              |                                     |
  ConditionalConstIteratorWithMask <------ ImageSpatialObject


Hope some of this is useful...

Regards -

  :: Gavin

-- 
Gavin Baker                                      Complex Systems Group
http://www.cs.mu.oz.au/~gavinb             The University of Melbourne