[Insight-developers] A fix for the bug

[Miller, James V (Research)]

I see Paul's point.  The GenerateInputRequestedRegion() for DiscreteGaussian
is asking for the proper amount of the original input.  The problem is that
the output requested region for the internal filter cannot set simple as the
output requested region.  So the GenerateInputRequestedRegion() is correct,
is the managing of the output requested region for the internal
mini-pipeline that is problematic.  Paul's solution of actually establishing
a complete mini-pipeline should work (I haven't looked at that code yet but
it should have the standard GraftOutput() calls and set the ReleaseData flag
on 2nd through Nth stage of the mini-pipeline.

Another solution is to manage the output requested region for the
mini-pipeline.  The original output of the DiscreteGaussian is grafted onto
the mini-pipeline.  This sets the region ivars for the mini-pipeline as well
as sets up the mini-pipeline to use the bulk data container of the
DiscreteGaussian output.  

What needs to be done is that each call to the mini-pipeline needs to have a
different output requested region.  The first call will be for the original
output requested region padded in N-1 of the dimensions (the mini-pipeline
will pad the last dimension for the input requested region).  The second
call will be for the original output requested region padded in N-2 of the
dimensions, etc. The last call for the original output requested region.

There are other places in the toolkit where we probably need to fix this.
Pretty much anywhere we are using a single filter in a loop where the filter
is NOT set up to require all the input and produce all the output.  

Setting up a full mini-pipeline as Paul suggested will also address this
problem. I think if the release data mechanism is used on the 2nd thru Nth
stage of the mini-pipeline, the two approaches will have identical memory
utilization.



-----Original Message-----
From: Paul Yushkevich [mailto:pauly at cognitica . com]
Sent: Wednesday, May 28, 2003 2:37 PM
To: insight-dev-list
Subject: [Insight-developers] A fix for the bug


Hi

I think I know the reason I was getting that bug (DiscreteGaussian 
filter crashing when using a single-slice RequestedRegion)
I think the problem is in the implementation of the GenerateData() 
method of the DG filter.  For a 3D image, this implementation applies 
NeighborhoodOperatorImageFilter 3 times.  First time it's applied to the 
input image, the second time it's applied to the output of the first 
attempt, the third time it's applied to the result of the second 
attempt. However, since there is only one filter object, I think (?) 
that the requested region is not being propagated from the output to the 
input correctly.  That is, the BufferedRegion of the output of the first 
application of the filter (call it BR1) ends up being equal to the 
RequestedRegion of the output of the entire DiscreteGaussian filter 
(call it RR).  That can't be right, because BR1 should include some 
padding, (i.e., BR1 = BR2 + padding, BR2 = BR3 + padding, BR3 = RR)

Here is a reimplementation of the DiscreteGaussian filter that 
constructs a whole mini-pipeline consisting of N filters with N 
operators.  This mini-pipeline ensures that the correct requested 
regions are propagated from the end of the pipeline to the input.  I am 
not sure if this introduces some inefficiencies into the code, and if 
there is a better way to fix the problem.

template< class TInputImage, class TOutputImage >
void
DiscreteGaussianImageFilter<TInputImage, TOutputImage>
::GenerateData()
{
  typename TOutputImage::Pointer output = this->GetOutput();
 
  output->SetBufferedRegion(output->GetRequestedRegion());
  output->Allocate();

  // Type definition for the internal neighborhood filter
  typedef NeighborhoodOperatorImageFilter<
    InputImageType, OutputImageType>              InternalFilterType;
  typedef typename InternalFilterType::Pointer InternalFilterPointer;

  // Create a chain of filters
  InternalFilterPointer *filter =
    new InternalFilterPointer[m_FilterDimensionality];

  // Create a series of operators
  GaussianOperator<OutputPixelType, ImageDimension> *oper =
    new 
GaussianOperator<OutputPixelType,ImageDimension>[m_FilterDimensionality];

  // Initialize and connect the filters
  for (unsigned int i = 0; i < m_FilterDimensionality; ++i)
    {
    // Set up the operator for this dimension
    oper[i].SetDirection(i);
    if (m_UseImageSpacing == true)
      {
      if (this->GetInput()->GetSpacing()[i] == 0.0)
        {
        itkExceptionMacro(<< "Pixel spacing cannot be zero");
        }
      else
        {
        oper[i].SetVariance(m_Variance[i] / 
this->GetInput()->GetSpacing()[i]);
        }
      }
    else
      {
      oper[i].SetVariance(m_Variance[i]);
      }

    oper[i].SetMaximumKernelWidth(m_MaximumKernelWidth);
    oper[i].SetMaximumError(m_MaximumError[i]);
    oper[i].CreateDirectional();

    // Set up the filter
    filter[i] = InternalFilterType::New();
    filter[i]->SetOperator(oper[i]);
    filter[i]->ReleaseDataFlagOn();

    // Connect the filter to the input or to the previous filter
    filter[i]->SetInput(i == 0 ? this->GetInput() : 
filter[i-1]->GetOutput());
    }

  // Graft this filters output onto the mini-pipeline so that the 
mini-pipeline
  // has the correct region ivars and will write to this filters bulk data
  // output.
  filter[m_FilterDimensionality-1]->GraftOutput( output );

  // Update the last filter in the chain
  filter[m_FilterDimensionality-1]->Update();
 
  // Graft the last output of the mini-pipeline onto this filters output so
  // the final output has the correct region ivars and a handle to the final
  // bulk data
  this->GraftOutput(output);

  // Clean up
  delete[] oper;
  delete[] filter;
}

-- 
--------------------------------
Paul A. Yushkevich, Ph.D.
President, Cognitica Corporation

17 Flemington Rd
Chapel Hill, NC 27517
Tel: 1-919-929-7652
--------------------------------


_______________________________________________
Insight-developers mailing list
Insight-developers at public . kitware . com
http://public . kitware . com/mailman/listinfo/insight-developers