[Insight-users] is there a problem in boundary face calculator ?
Kristen Zygmunt
krismz at sci.utah.edu
Wed May 11 17:25:15 EDT 2011
Hello,
I was noticing a problem with the boundary face calculator where some pixels were visited multiple times.
I think this could be a real problem, especially for multithreading where it is assumed that images are divided up into non-overlapping regions. I did some digging in the threads and came across two discussions, one reported the problem on October 17, 2010 (first conversation below) and another provided a fix on November 15, 2010 (second conversation below). I have taken the fix and applied it to the 3.20 branch, that patch is attached. Also, I have written a simple test driver that demonstrates the problem and attached it as well.
In response to the first conversation below, this is a problem for all images. Specifically, my test driver creates a 3D image that defaults to size 4x4x4, though the size can be passed in via an argument and uses a neighborhood radius of 1. I have used the test to directly confirm this is an issue for 3x3x3, 4x4x4, and 10x10x10 images.
Should I create a bug report for this, or simply create a gerrit patch for ITK4?
Thanks!
Kris
Conversation reporting the problem:
> Date: Sun, 17 Oct 2010 08:16:14 -0400
> From: Luis Ibanez <luis.ibanez at kitware.com>
> Subject: Re: [Insight-users] is there a problem in boundary face
> calculator ?
> To: devieill at irit.fr
> Cc: insight-users at itk.org
> Message-ID:
> <AANLkTimnydwJcaoaUL=PGJVDmCw0++hyckdPyCdOJFGb at mail.gmail.com>
> Content-Type: text/plain; charset="utf-8"
>
> Hi Francois,
>
> Thanks for pointing this out.
>
> You may have found a bug indeed.
>
> The assumption of the boundary face
> calculator is that the radius of the
> neighborhood is smaller than half
> of the image size.
>
> Let me attempt to illustrate this in a
> diagram:
>
> XXXXXXXXXXXXXXXXXXXX
> XXXXXXXXXXXXXXXXXXXX
> HHHOOOOOOOOOOOOOOHHH
> HHHOOOOOOOOOOOOOOHHH
> HHHOOOOOOOOOOOOOOHHH
> HHHOOOOOOOOOOOOOOHHH
> XXXXXXXXXXXXXXXXXXXX
> XXXXXXXXXXXXXXXXXXXX
>
> Here we have the top and
> bottom faces for a neighborhood
> whose radius is (3,2) (3 in the
> horizontal direction, and 2 in the
> vertical direction).
>
> Therefore, the "H" colums
> form faces of width=3
>
> while the "X" rows form faces
> of height=2
>
> Leaving ample room for the
> central "face" marked by the
> "O" symbols.
>
> If we start increasing the values
> of the radius, the "H" and "X"
> regions will increase in size,
> at the same time that the "O"
> (the inside face) will decrease
> in size.
>
> Therefore, it is possible to get
> to a point where the "O" region
> is empty, or even further, where
> the left and right "H" regions
> start overlapping.
>
> In your particular case,
>
> A) You have an image of 4x4 pixels.
>
> B) You set the radius to 1,1
> therefore to a neighborhood
> size of 3x3
>
> C) You ask the calculator to partition
> a subregion of the image, defined
> by
>
> Index = (1,0)
> Size = (1,2)
>
>
> So, in practice, you "image of interest''
> is only of size 1x2.
>
>
> Should you have used the full size of
> the 4x4 image, you should have found
> faces as in the diagram below:
>
>
> BBBB
> CAAD
> CAAD
> EEEE
>
> Where "A" is the inside face,
> and "B","C","D","E", are the lateral
> faces
>
> Since in (3) you limited the region
> to only the "O"s in the diagram below:
>
> XOOX
> XXXX
> XXXX
> XXXX
>
> Then the partition of the "OO" region
> into faces becomes quite challenging.
>
> The resulting "inside" region is certainly
> correct in retuning an empty size.
>
> As for the other regions, I agree with
> you in that they are not consistent...
>
> However,
> I'm having trouble figuring out what
> the "correct" answer should have been...
>
>
> I would be tempted to say that the
> request is not solvable, and that the
> calculator should have simply thrown
> and exception here.
>
>
> Would you consider that a reasonable
> behavior for the face calculator ?
>
>
>
> Regards,
>
>
> Luis
>
>
> -----------------------------------------------------
> On Fri, Oct 1, 2010 at 12:20 PM, <devieill at irit.fr> wrote:
>
>> Hello all,
>> the itk documentation state that the boundary face calculator should output
>> 2N+1 region where N is considered to be the region/image dimension.
>> However when creating a itk::Image<double,2> containing 4 pixels and
>> running the following simple code :
>> ____________________________________
>>
>> typedef itk::NeighborhoodAlgorithm::ImageBoundaryFacesCalculator<
>> Image2D > FaceCalculatorType;
>> itk::NeighborhoodAlgorithm::ImageBoundaryFacesCalculator< Image2D>
>> faceCalculator;
>> FaceCalculatorType::FaceListType faceList;
>> itk::NeighborhoodAlgorithm::ImageBoundaryFacesCalculator<
>> Image2D>::FaceListType::iterator fit;
>>
>> Image2D::RegionType region;
>> Image2D::IndexType index;
>> Image2D::SizeType size;
>> itk::Size<2> radius;
>> index[0] = 1;
>> index[1] = 0;
>> size[0] = 1;
>> size[1] = 2;
>> region.SetIndex(index);
>> region.SetSize(size);
>> radius[0]=1;
>> radius[1]=1;
>> faceList = faceCalculator(img1, region, radius);
>> std::cout << " face list number " << faceList.size() << std::endl ;
>> for (fit = faceList.begin(); fit!= faceList.end(); ++fit) {
>> std::cout << " current face is " << *fit << std::endl ;
>> }
>> _______________________________
>>
>>
>> I have the following output :
>>
>> ______________________________
>> face list number 4
>> current face is ImageRegion (0x25c0f10)
>> Dimension: 2
>> Index: [1, 1]
>> Size: [0, 0]
>>
>> current face is ImageRegion (0x25c0f50)
>> Dimension: 2
>> Index: [1, 0]
>> Size: [1, 2]
>>
>> current face is ImageRegion (0x25c0f90)
>> Dimension: 2
>> Index: [1, 0]
>> Size: [1, 1]
>>
>> current face is ImageRegion (0x25c0fd0)
>> Dimension: 2
>> Index: [1, 1]
>> Size: [1, 1]
>> ____________________________________
>>
>>
>> Despite the fact that there is NOT the right number of regions,
>> what bothers me the most is that one gets duplicated pixels, or to put it
>> simply the regions given by the algorithm are overlapping.
>> Assuming that one build a filter using neighboorhood filters this
>> become troublesome (imagine with multithreading...).
>> My question are :
>> 1) is the boundary calculator broken ?
>> 2) what should be used instead to design fast filters requiring a radius ?
>> 3) Can one have neuman boundary conditions with it ?
>>
>> Regards,
>> de Vieilleville Fran?ois
>>
>>
>
>
Conversation fixing the problem:
> Date: Mon, 15 Nov 2010 17:41:28 +1100
> From: wanlin <wanlinzhu at gmail.com>
> Subject: Re: [Insight-users] ShapedNeighborhoodIterator problem
> To: Dawood Masslawi <masslawi at gmail.com>
> Cc: insight-users at itk.org
> Message-ID:
> <AANLkTinshsf2reaCTjSXHRh61dAp1OVXW0aWq9493jyN at mail.gmail.com>
> Content-Type: text/plain; charset="iso-8859-1"
>
> For issue 5). I ran cross it before and make some changes in
> itkNeighborhoodAlgorithm.txx. Here is the patch, For 5x5 image with kernel
> size 3x3, the output of regions list are
>
> 9
> ImageRegion (0xb8b830)
> Dimension: 2
> Index: [1, 1]
> Size: [3, 3]
>
> 5
> ImageRegion (0xb8b870)
> Dimension: 2
> Index: [0, 0]
> Size: [1, 5]
>
> 5
> ImageRegion (0xb8b8b0)
> Dimension: 2
> Index: [4, 0]
> Size: [1, 5]
>
> 3
> ImageRegion (0xb8b8f0)
> Dimension: 2
> Index: [1, 0]
> Size: [3, 1]
>
> 3
> ImageRegion (0xb8b930)
> Dimension: 2
> Index: [1, 4]
> Size: [3, 1]
>
>
>
>
> On Mon, Nov 15, 2010 at 11:37 AM, Dawood Masslawi <masslawi at gmail.com>wrote:
>
>>
>> - 1) If I'm going to do this type of thing:
>> -
>> - IteratorType::OffsetType top = {{0,-1}};
>> - IteratorType::OffsetType bottom = {{0,1}};
>> - IteratorType::OffsetType left = {{-1,0}};
>> - IteratorType::OffsetType right = {{1,0}};
>> -
>> - Then get the pixels with:
>> - iterator[top][0]
>> -
>> - is there an advantage to using the ShapedNeighborhoodIterator over
>> just a NeighborhoodIterator?
>>
>> ------------------------------
>>
>> You can use offsets for both ShapedNeighborhoodIterator and regular
>>
>> NeighborhoodIterator, the difference is in performance. A regular
>>
>> NeighborhoodIterator would need more pixels to form a rectangular
>>
>> neighborhood (which artificially could be considered as a
>>
>> "ShapedNeighborhoodIterator" using offsets) but a
>>
>> ShapedNeighborhoodIterator would require less pixels in the neighborhood
>>
>> to form a shaped neighborhood, less pixels mean less memory allocation
>>
>> and less boundary checking which could remarkably improve the performance
>>
>> specially for large neighborhoods in large and multidimensional images
>> plus,
>>
>> handling shaped neighborhoods is easier with
>> the ShapedNeighborhoodIterator.
>>
>> ------------------------------
>>
>>
>> - 2) If I do want to use the ShapedNeighborhood style, is this correct:
>> -
>> - IteratorType::IndexListType::const_iterator indexIterator =
>> iterator.GetActiveIndexList().begin();
>> - while (indexIterator != iterator.GetActiveIndexList().end())
>> - {
>> - std::cout << (int)iterator[*indexIterator][0] << " ";
>> - ++indexIterator;
>> - }
>>
>>
>> ------------------------------
>> Yes, the following is also a good practice:
>>
>> IteratorType::ConstIterator ci;
>>
>> for (ci = iterator.Begin(); ci != iterator.End(); ci++)
>>
>> {
>>
>> std::cout << ci.get() << std::endl;
>>
>> }
>>
>> ------------------------------
>>
>>
>> - 3) If I loop over the ActiveIndexList using the FaceCalculator, is
>> the idea that in all regions except the 0th region I need to do something
>> like:
>> -
>> - IteratorType::IndexListType::const_iterator indexIterator =
>> iterator.GetActiveIndexList().begin();
>> - while (indexIterator != iterator.GetActiveIndexList().end())
>> - {
>> - bool valid;
>> - iterator.GetPixel(*indexIterator,valid);
>> - if(valid) { do something }
>> - ++indexIterator;
>> - }
>> -
>> - Alternatively, I could use 4 separate loops (one for each face
>> region) where I know which neighbor isn't valid. I would have to modify the
>> ActiveIndexList to reflect this missing pixel before each loop.
>>
>> ------------------------------
>>
>> Usually for neither of face calculator or neighborhood iterator you need to
>> perform
>>
>> the boundary check yourself (notice that the iterator itself has to iterate
>> over the
>>
>> boundary to check for its validity), the face calculator splits the image
>> into boundary
>>
>> and non-boundary regions for the iterator to not to check for all of the
>> pixels and
>>
>> only perform the boundary checking for pixels in the face region (pixels
>> with distance
>>
>> from the image boundary equal or less than the neighborhood radius). This
>> would also
>>
>> improve the performance specially for large datasets.
>>
>> ------------------------------
>>
>> - 4) I believe
>> - iterator[top][0]
>> - is equivalent to but faster than iterator.GetPixel(top). Is this
>> correct?
>>
>> ------------------------------
>> I'm not sure about this since I haven't seen a significant performance
>> difference.
>>
>>
>> -
>> ------------------------------
>> - 5) Am I correct that the face calculator returns overlapping regions?
>> If I use a 5x5 image with a 3x3 kernel, it returns regions of size
>> - 9
>> - 5
>> - 5
>> - 5
>> - 5
>> -
>> - totaling 29, when there are only 25 pixels. I don't see any functions
>> to turn on/off this overlap. How is this usually handled?
>>
>> ------------------------------
>>
>> To the best of my knowledge this is still an unsolved issue which was
>> brought up
>>
>> by other users in the list as well. Usually the face calculator assumes
>> that the
>>
>> face region is far smaller than the non-boundary region which could be an
>> issue
>>
>> for small images or very large neighborhoods.
>>
>>
>> Regards,
>>
>> Dawood
>>
>>
>> _____________________________________
>> Powered by www.kitware.com
>>
>> Visit other Kitware open-source projects at
>> http://www.kitware.com/opensource/opensource.html
>>
>> Kitware offers ITK Training Courses, for more information visit:
>> http://www.kitware.com/products/protraining.html
>>
>> Please keep messages on-topic and check the ITK FAQ at:
>> http://www.itk.org/Wiki/ITK_FAQ
>>
>> Follow this link to subscribe/unsubscribe:
>> http://www.itk.org/mailman/listinfo/insight-users
>>
>>
> -------------- next part --------------
> An HTML attachment was scrubbed...
> URL: <http://www.itk.org/pipermail/insight-users/attachments/20101115/93d58843/attachment-0001.htm>
> -------------- next part --------------
> A non-text attachment was scrubbed...
> Name: itkNeighborhoodAlgorithm.txx.patch
> Type: text/x-patch
> Size: 2375 bytes
> Desc: not available
> URL: <http://www.itk.org/pipermail/insight-users/attachments/20101115/93d58843/attachment-0001.bin>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20110511/8f793867/attachment.htm>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: ITK_3.20.neighborhood.patch
Type: application/octet-stream
Size: 2762 bytes
Desc: not available
URL: <http://www.itk.org/pipermail/insight-users/attachments/20110511/8f793867/attachment.obj>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20110511/8f793867/attachment-0001.htm>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: itkNeighborhoodAlgorithmTest.cxx
Type: application/octet-stream
Size: 4559 bytes
Desc: not available
URL: <http://www.itk.org/pipermail/insight-users/attachments/20110511/8f793867/attachment-0001.obj>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20110511/8f793867/attachment-0002.htm>
More information about the Insight-users
mailing list