[Insight-developers] Re: [Insight-users] One error and one is sue with itkImagePCAShapeModelEstimator

[Miller, James V (Research)]

Zachary,

The RequestedRegion is what filters need to operate on. The pipeline
(and filters) are designed so that the RequestedRegion should always
be a subset of the BufferedRegion. Since your GenerateInputRequestedRegion
sets the RequestedRegion to the LargestPossibleRegion, it 
should be safe to use the RequestedRegion (and not the BufferedRegion
in the CalculateInnerProducts).

What your estimator needs is for the calculation of the PCA to occur
over regions of the inputs that are the same size.  This means the
GenerateInputRequestedRegion method needs to guarentee the inputs have
the same RequestedRegion size.

Currently, you are driving this from the "first" input.  This is a 
reasonable choice.  Another option is to scan all the inputs to find the 
smallest input and make that the master. Then it becomes a question of
how to how to change the start indices if you end up changing the size.

The cleanest thing to do, however, is to probably drive this off an 
output requested region. The outputs should be images so this seems
reasonable.  The input RequestedRegions are then just copies of "an"
output RequestedRegion. I think the default pipeline implementations 
should handle this appropriately (meaning the filter probably does not
need a GenerateInputRequestedRegion method or an
EnlargeOutputRequestedRegion
method).  The PCA decomposition will be on whatever RequestedRegion 
is asked for on output.  If someone calls UpdateLargestPossibleRegion()
on the estimator, then the PCA decomposition will be over the entire image.

This will displace the complexity to the user.  Where they will use a 
series of ExtractImageFilters, RegionOfInterestImageFilters, and 
ChangeInformationFilters to get all the inputs a consistent size... If 
we do not displace the complexity to the user, then the filter needs
to know which portions of each input to use in the PCA.  For instance, 
lets say I have a number of images, each having a small object that 
I want to model.  In medical image analysis, it is rare that the object 
of interest will appear at the same spot in each  image.  Therefore, I would

expect to have to align and crop a block around the object to send it to the
PCA.
Adding this complexity to the PCA filter seems unnecessary.  The PCA filter
should probably just check the inputs are the same size.


Jim


-----Original Message-----
From: Zachary Pincus [mailto:zpincus at stanford.edu]
Sent: Sunday, May 16, 2004 7:02 AM
To: Luis Ibanez; insight-developers at itk.org
Subject: [Insight-developers] Re: [Insight-users] One error and one
issue with itkImagePCAShapeModelEstimator


> 1) It seems reasonable to allow the ImagePCAShapeModelEstimator
>    to accept images of different origin, different StartIndex
>    and only enforce the pixel spacing and the image size (in
>    pixel along every dimension) to be the same among all the
>    input images.
>
>    Please send us a patch...

I've included a new implementation of the GenerateInputRequestedRegion 
method of itkPCAShapeModelEstimatorClass so images with regions with 
different indices, but the same size (and spacing) can work.

Now, the original behavior (which I did not change) is to set all input 
requested regions to the largest possible region, and claim to use 
that. However, elsewhere in the code (CalculateInnerProducts), the 
buffered region of the inputs is what is actually used...

Perhaps it would be better if this were all done with the requested 
region? If the user bothered to set them, why blow them away? Or would 
that break too much code that relies on not having to set a requested 
region?

Anhyow, the diff is pretty messy since this is basically a different 
function, so I've just pasted in the new function.

Also attached is a modified test case that exercises this new 
functionality, and diffs to the old test.


Zach Pincus

Department of Biochemistry and Program in Biomedical Informatics
Stanford University School of Medicine


New GenerateInputRequestedRegion:

/**
  * Requires all of the inputs to have a LargestPosibleRegion and 
Spacing the of
  * same size.
  */
template<class TInputImage, class TOutputImage>
void
ImagePCAShapeModelEstimator<TInputImage,TOutputImage>
::GenerateInputRequestedRegion()
{
   Superclass::GenerateInputRequestedRegion();
   // Calling this causes problems because it sets the input requested 
regions
   // to various output requested regions, which might not overlap, 
since the inputs
   // and outputs aren't very closely linked in the PCA case.
   // If we REALLY want to call this, then some heroics are necessary to 
undo
   // the damage it does. Specifically, we will have to set all the 
inputs'
   // requested regions back to *their* LargestPossibleRegion, and not 
that of
   // some random output.
	
   if ( this->GetInput(0) )
     {

     // Get the requested region and spacing of the first input
	this->GetInput(0)->SetRequestedRegionToLargestPossibleRegion();
	typename TInputImage::RegionType masterRequestedRegion =
	  this->GetInput(0)->GetRequestedRegion();
     typename TInputImage::SpacingType masterSpacing =
       this->GetInput(0)->GetSpacing();

     // Set the requested region of the remaining input to their largest 
possible
     // regions, and make sure the sizes and spacings are proper.
     unsigned int idx;
     for (idx = 1; idx < this->GetNumberOfInputs(); ++idx)
       {
       if ( this->GetInput(idx) )
         {
		
         
this->GetInput(idx)->SetRequestedRegionToLargestPossibleRegion();
		
		typename TInputImage::RegionType requestedRegion =
           this->GetInput(idx)->GetRequestedRegion();
         typename TInputImage::SpacingType spacing =
		  this->GetInput(0)->GetSpacing();
         if( masterRequestedRegion.GetSize() != 
requestedRegion.GetSize() ||
			masterSpacing != spacing)
           {
           itkExceptionMacro("Size and spacing of input " << idx <<
             " is not the same as that of input 0" );
           }

         } // if ( this->GetIntput(idx))
       } // for idx
     } // if( this->GetInput(0) )

}


Diffs for new and old test code 
(itkImagePCAShapeModelEstimatorTest.cxx):
< Old
 > New
83c83
<   InputImageType::IndexType index;
---
 >   InputImageType::IndexType index, index2;
85c85,87
<   InputImageType::RegionType region;
---
 >   index2.Fill(4);
 >
 >   InputImageType::RegionType region, region2;
89a92,94
 >   region2.SetSize( inputImageSize );
 >   region2.SetIndex( index2 );
 >
105,106c110,111
<   image2->SetLargestPossibleRegion( region );
<   image2->SetBufferedRegion( region );
---
 >   image2->SetLargestPossibleRegion( region2 );
 >   image2->SetBufferedRegion( region2 );