[Insight-users] Pipeline question
Paul Yushkevich
pauly at cognitica . com
Thu, 15 May 2003 10:55:37 -0400
This is a multi-part message in MIME format.
--------------040307070300060601000703
Content-Type: text/plain; charset=ISO-8859-1; format=flowed
Content-Transfer-Encoding: 7bit
Miller, James V (Research) wrote:
I wrote my own filter to accomodate the needs of our application, and
here is the code for it. I've marked the place of interest with a large
//**** comment.
Paul.
>Paul,
>
>What filter are you using to extract the slices?
>
>Jim
>
>
>-----Original Message-----
>From: Paul Yushkevich [mailto:pauly@cognitica.com]
>Sent: Thursday, May 15, 2003 9:35 AM
>To: Miller, James V (Research)
>Cc: insight
>Subject: Re: [Insight-users] Pipeline question
>
>
>James,
>
>Let me try to clarify the problem.
>
>Let's say I have a large 3D image that never changes. I then have the 3
>slicing filters, each of which is followed by an intensity remapping
>filter, then followed by a padding filter that makes the slices 256x256
>and then by another 'mapper' that stores the slices in texture memory
>for display. I update these mappers each time that my GUI needs to
>redraw. However, I would like the pipeline to cache the output of the
>slicing and preprocessing and texturing filters as long as I don't
>change their parameters (e.g. the slice number to display)
>
>What happens is that this display pipeline gets executed every time that
>Update is called on the texturing mappers at the end of the pipeline. I
>think that that happens because the call to Update causes the
>RequestedRegion to be propagated backwards through the pipeline. For
>each of the 3 slices, a different RequestedRegion is set in the 3D
>image. As a result, the 3D Image's MTime gets incremented for each of
>these Updates(), and then the whole pipeline gets executed, and no
>caching occurs.
>
>Now, I can get around this problem by having the slicing filters set
>their input's requested region to be the LargestPossibleRegion. That's
>fine for a static 3D image that I have described, but imagine that the
>3D image is an output of another pipeline, and I want to run that
>pipeline only for the pixels in the slices. In this situation, setting
>RequestedRegion to the LargestPossibleRegion is not what I want
>
>Thank you!
>
>Paul.
>
>
>Miller, James V (Research) wrote:
>
>
>
>>Paul,
>>
>>Let me see if I understand the problem correctly.
>>
>>You have a pipeline where near the end of the pipeline you have 3 separate
>>slice extraction filters. When you call Update() on the slice extraction
>>filters in succession, the entire pipeline re-executes? Or just the slice
>>extraction filters re-execute?
>>
>>Is the input image to the three slice extraction filters "big enough" to
>>support all three requests? If every filter before the slice extraction
>>filter that you are executing properly handles streaming, then filter
>>
>>
>before
>
>
>>the slice extraction filters would probably have only produced just enough
>>output to support the first slice extraction request. When the second
>>request arrives, the pipeline would re-execute to generate the new
>>
>>
>requested
>
>
>>region.
>>
>>One thing you might want to check is whether the BufferedRegion is changing
>>on the input to the slice extraction filters after each successive call to
>>Update() on the slice extraction filters. If it is changing, then it means
>>that when the filter before the slice extraction gets each request, it
>>thinks it needs to re-execute. One thing you can try is to call
>>UpdateLargestPossibleRegion() on the filter before the slice extraction
>>filters so that all the output (to that point) is available. Then the
>>
>>
>slice
>
>
>>extraction filters should never cause the filter before the slice
>>
>>
>extraction
>
>
>>filter to re-execute.
>>
>>The last time I tested things like this explictly was a couple of months
>>
>>
>ago
>
>
>>when I added InPlace filters. Here I exercised the pipeline quite a bit to
>>make sure branches of a pipeline properly executed. But your issue may be
>>particular to the slice extraction filter.
>>
>>But if the filter before the slice extraction filter is not re-executing
>>
>>
>but
>
>
>>the slice extraction filters are always re-executing, then that is a
>>different issue. This could therefore be an MTime issue.
>>
>>If you could whittle this down to a small test program that would certainly
>>help us out.
>>
>>
>>
>>-----Original Message-----
>>From: Paul Yushkevich [mailto:pauly@cognitica.com]
>>Sent: Wednesday, May 14, 2003 5:02 PM
>>To: insight
>>Subject: [Insight-users] Pipeline question
>>
>>
>>Hi,
>>
>>I have a small question about the pipeline.
>>
>>I have a large image that is the input to three slice extraction filters
>>that operate on different regions of the image. Now, every time I call
>>Update() on the outputs of these filters, the entire pipeline gets
>>rerun, even though neither the input image nor the parameters of the
>>filters had changed.
>>
>>I traced this down to the fact that the input image is Modified() every
>>time its requested region changes. As different RequestedRegions are
>>propagated from the three different filters, the input image's MTime
>>gets incremented, and the image appears as if it has changed, even
>>though it really hadn't.
>>
>>Is there any way around this issue? I want to be able to call Update()
>>frequently on the outputs of the three slice filters, but only have the
>>filters run if the image or the filter parameters have changed.
>>
>>Paul
>>
>>
>>
>>
>>
>
>
>
>
--
--------------------------------
Paul A. Yushkevich, Ph.D.
President, Cognitica Corporation
17 Flemington Rd
Chapel Hill, NC 27517
Tel: 1-919-929-7652
--------------------------------
--------------040307070300060601000703
Content-Type: text/plain;
name="IRISSlicer.h"
Content-Transfer-Encoding: 7bit
Content-Disposition: inline;
filename="IRISSlicer.h"
#ifndef __IRISSlicer_h_
#define __IRISSlicer_h_
#include <ImageCoordinateTransform.h>
#include <itkImageToImageFilter.h>
#include <itkImageSliceConstIteratorWithIndex.h>
#include <itkImageRegionIterator.h>
/**
* A slice extraction filter for 3D images. This filter takes a transform
* from image space (x=pixel, y=line, z=slice) to display space. A slice
* shows either x-y, y-z, or x-z in the display space. This filter is necessary
* because the origin in the slice can correspond to any corner of the image, not
* just to the origin of the image. This filter can traverse the image in different
* directions, accomodating different positions of the display space origin in the
* image space.
*/
template <class TPixel>
class ITK_EXPORT IRISSlicer
: public itk::ImageToImageFilter<itk::Image<TPixel,3>, itk::Image<TPixel,2> >
{
public:
/** Standard class typedefs. */
typedef IRISSlicer Self;
typedef itk::Image<TPixel,3> InputImageType;
typedef itk::Image<TPixel,2> OutputImageType;
typedef itk::ImageToImageFilter<InputImageType,OutputImageType> Superclass;
typedef itk::SmartPointer<Self> Pointer;
typedef itk::SmartPointer<const Self> ConstPointer;
/** Method for creation through the object factory. */
itkNewMacro(Self);
/** Run-time type information (and related methods). */
itkTypeMacro(IRISSlicer, ImageToImageFilter);
/** Some typedefs. */
typedef typename InputImageType::ConstPointer InputImagePointer;
typedef typename InputImageType::RegionType InputImageRegionType;
typedef typename OutputImageType::Pointer OutputImagePointer;
typedef typename OutputImageType::RegionType OutputImageRegionType;
typedef itk::ImageSliceConstIteratorWithIndex<InputImageType> InputIteratorType;
typedef itk::ImageRegionIterator<OutputImageType> OutputIteratorType;
/** Set the transform that maps the image space into the display space */
void SetImageToDisplayTransform(const ImageCoordinateTransform &imageToDisplay);
/** Get the transforms stored in the slicer */
itkGetMacro(ImageToDisplay,ImageCoordinateTransform);
itkGetMacro(DisplayToImage,ImageCoordinateTransform);
/** Set the current slice index */
itkSetMacro(SliceIndex,unsigned int);
itkGetMacro(SliceIndex,unsigned int);
/** Set the current axis direction */
itkSetMacro(SliceAxis,unsigned int);
itkGetMacro(SliceAxis,unsigned int);
protected:
IRISSlicer();
virtual ~IRISSlicer() {};
/**
* IRISSlicer can produce an image which is a different
* resolution than its input image. As such, IRISSlicer
* needs to provide an implementation for
* GenerateOutputInformation() in order to inform the pipeline
* execution model. The original documentation of this method is
* below.
*
* \sa ProcessObject::GenerateOutputInformaton() */
virtual void GenerateOutputInformation();
/**
* This method maps an input region to an output region
*/
virtual void CallCopyOutputRegionToInputRegion(InputImageRegionType &destRegion,
const OutputImageRegionType &srcRegion);
/**
* IRISSlicer is not implemented as a multithreaded filter.
* \sa ImageToImageFilter::GenerateData()
*/
virtual void GenerateData();
private:
IRISSlicer(const Self&); //purposely not implemented
void operator=(const Self&); //purposely not implemented
// Transform from image to the slices
ImageCoordinateTransform m_ImageToDisplay;
// Transform from slices to the image
ImageCoordinateTransform m_DisplayToImage;
// Current slice in each of the dimensions
unsigned int m_SliceIndex;
// Current axis direction
unsigned int m_SliceAxis;
// Axes in display space that correspond to the z,x,y axes of the slice
// (The z axes of the slice is the dimension in which the slice moves)
Vector3i m_DisplayAxes;
// Axes in image space that correspond to the z,x,y axes of the slice
// (The z axes of the slice is the dimension in which the slice moves)
Vector3i m_ImageAxes;
// The direction in which the image axes are traversed (+1 or -1)
Vector3i m_AxesDirection;
// This method computes the axes. Since this is so fast, the method is run
// each time the axes are needed
void ComputeAxes();
// The worker method in this filter
void CopySlice(InputIteratorType itImage, OutputIteratorType itSlice,
int sliceDir, int lineDir);
};
template<class TPixel>
IRISSlicer<TPixel>
::IRISSlicer()
{
// There is a single input to the filter
SetNumberOfRequiredInputs(1);
// There are three outputs from the filter
SetNumberOfRequiredOutputs(1);
// Initialize slice indices
m_SliceIndex = 0;
m_SliceAxis = 0;
}
template<class TPixel>
void IRISSlicer<TPixel>
::GenerateOutputInformation()
{
// Get pointers to the inputs and outputs
typename Superclass::InputImageConstPointer inputPtr = GetInput();
typename Superclass::OutputImagePointer outputPtr = GetOutput();
// The inputs and outputs should exist
if (!outputPtr || !inputPtr) return;
// Get the size and origin of the output image
InputImageRegionType inputRegion = inputPtr->GetLargestPossibleRegion();
itk::Size<3> inputSize = inputRegion.GetSize();
Vector3d inputSpacing(inputPtr->GetSpacing());
// Compute the axes directions
ComputeAxes();
// Create the new size for the region
itk::Size<2> sliceSize;
sliceSize[0] = inputSize[m_ImageAxes(1)];
sliceSize[1] = inputSize[m_ImageAxes(2)];
// Set the region of the output slice
outputPtr->SetRegions(sliceSize);
// Compute the origin of the new image
Vector2d sliceOrigin(0.0);
outputPtr->SetOrigin(sliceOrigin.data_block());
// Compute the spacing of the new image
Vector2d sliceSpacing(inputSpacing(m_ImageAxes(1)),
inputSpacing(m_ImageAxes(2)));
outputPtr->SetSpacing(sliceSpacing.data_block());
}
template<class TPixel>
void IRISSlicer<TPixel>
::CallCopyOutputRegionToInputRegion(InputImageRegionType &destRegion,
const OutputImageRegionType &srcRegion)
{
// Start with the largest possible region
destRegion = GetInput()->GetLargestPossibleRegion();
// *************************************************
// This is the code I had to comment out because it
// causes the pipeline to be re-evaluated each time
// *************************************************
/*
// Compute the axes directions
ComputeAxes();
// Set the size of the region to 1 in the requested direction
destRegion.SetSize(m_ImageAxes(0),1);
// Set the index of the region in that dimension to the number of the slice
destRegion.SetIndex(m_ImageAxes(0),m_SliceIndex);
*/
}
template<class TPixel>
void IRISSlicer<TPixel>
::GenerateData()
{
// Here's the input and output
InputImagePointer inputPtr = GetInput();
OutputImagePointer outputPtr = GetOutput();
// Do we really need this?
// ProgressReporter progress(this, threadId, outputRegionForThread.GetNumberOfPixels());
// Allocate (why is this necessary?)
outputPtr->Allocate();
// Compute the region in the image for which the slice is being extracted
InputImageRegionType inputRegion = inputPtr->GetLargestPossibleRegion();
// Set the size of the region to 1 in the requested direction
inputRegion.SetSize(m_ImageAxes(0),1);
// Set the index of the region in that dimension to the number of the slice
inputRegion.SetIndex(m_ImageAxes(0),m_SliceIndex);
// Compute the axes directions
ComputeAxes();
// Create an iterator that will parse the desired slice in the image
InputIteratorType it(inputPtr,inputRegion);
it.SetFirstDirection(m_ImageAxes(1));
it.SetSecondDirection(m_ImageAxes(2));
// Create an iterator into the slice itself
OutputIteratorType itSlice(outputPtr,outputPtr->GetLargestPossibleRegion());
// Copy the contents
CopySlice(it,itSlice,m_AxesDirection(1),m_AxesDirection(2));
cout << "SLICING " << m_SliceAxis << "; Type " << sizeof(TPixel) << endl;
}
template<class TPixel>
void IRISSlicer<TPixel>
::ComputeAxes()
{
// This is how the axes in the display space map onto the U,V axes of the
// individual slices
static unsigned int linkedAxes[3][3] = {{0,1,2},{1,0,2},{2,0,1}};
// Compute the slice axes in display and image
for(unsigned int i=0;i<3;i++)
{
// Compute the directions of the slice in display space
m_DisplayAxes(i) = linkedAxes[m_SliceAxis][i];
// Map the dimensions into the image space
m_ImageAxes(i) =
m_DisplayToImage.GetCoordinateIndexZeroBased(m_DisplayAxes(i));
// Compute the direction of traversal in these dimensions
m_AxesDirection(i) =
m_DisplayToImage.GetCoordinateOrientation(m_DisplayAxes(i));
}
}
template<class TPixel>
void IRISSlicer<TPixel>
::SetImageToDisplayTransform(const ImageCoordinateTransform &imageToDisplay)
{
// Set the transform
m_ImageToDisplay = imageToDisplay;
m_DisplayToImage = imageToDisplay.Inverse();
// Modified!
Modified();
}
template<class TPixel>
void IRISSlicer<TPixel>
::CopySlice(InputIteratorType itImage, OutputIteratorType itSlice,
int sliceDir, int lineDir)
{
if (sliceDir > 0 && lineDir > 0) {
itImage.GoToBegin();
while (!itImage.IsAtEndOfSlice()) {
while (!itImage.IsAtEndOfLine()) {
itSlice.Set(itImage.Value());
++itSlice;
++itImage;
}
itImage.NextLine();
}
} else if (sliceDir > 0 && lineDir < 0) {
itImage.GoToBegin();
while (!itImage.IsAtEndOfSlice()) {
itImage.NextLine();
itImage.PreviousLine();
while (1) {
itSlice.Set(itImage.Value());
++itSlice;
if (itImage.IsAtReverseEndOfLine())
break;
else
--itImage;
}
itImage.NextLine();
}
} else if (sliceDir < 0 && lineDir > 0) {
itImage.GoToReverseBegin();
while (!itImage.IsAtReverseEndOfSlice()) {
itImage.PreviousLine();
itImage.NextLine();
while (!itImage.IsAtEndOfLine()) {
itSlice.Set(itImage.Value());
++itSlice;
++itImage;
}
itImage.PreviousLine();
}
} else if (sliceDir < 0 && lineDir < 0) {
itImage.GoToReverseBegin();
while (1) {
while (1) {
itSlice.Set(itImage.Value());
++itSlice;
if (itImage.IsAtReverseEndOfLine())
break;
else
--itImage;
}
if (itImage.IsAtReverseEndOfSlice())
break;
else
itImage.PreviousLine();
}
}
}
#endif
--------------040307070300060601000703--