Hough Transform 2D Circles Image Filter GetCircles patch.


(Niels Dekker) #21

Thanks for your reply, Matt!

So filter->GetCircles() only actively searches for circles when the
’MTime’ of the filter has changed. (Otherwise it just returns the
previously found circles.)

When an input image of a filter is modified (by image->Modified()),
and filter->Update() is called, do you think that in general,
filter->GenerateData() should take care of updating the MTime?

Should HoughTransform2DCirclesImageFilter::GenerateData() be fixed to
call this->Modified()?


(Matt McCormick) #22

Good question.

->Modified() does not need to be called in a filter’s GenerateData() method because the pipeline system takes care of this.


(Niels Dekker) #23

Thanks, Matt. However… with the following code, I expected each filter->GetCircles() to do a new search for circles, on each iteration:

filter->SetInput(image);

for (int i = 0; i < n; ++i)
{
    image->FillBuffer(n)
    image->Modified();
    filter->Update();
    filter->GetCircles();
}

Unfortunately that did not happen, because filter->GetCircles() only searches for new circles when this->GetMTime() has changed. Otherwise it just returns the previously found circles:

Do you think in this case the user code is wrong?


(Matt McCormick) #24

Perhaps m_OldModifiedTime should be updated in GenerateData().


(Bradley Lowekamp) #25

Proper execution based in on modified times is a tricky thing. It it best handled in the pipeline.

What ever original purpose GetCircles had for doing computation outside of the pipeline line update is no longer there. Perhaps it was related to the former n argument.

I think the computation done in the GetCircle method should be move to be part of the pipeline update. Likely in the AfterThreadedGenerateData.


(Matt McCormick) #26

That’s a good idea.


(Dženan Zukić) #27

Will you do it, @Niels_Dekker?


(Niels Dekker) #28

If you would like to do it @dzenanz please go ahead! :slight_smile:

Honestly, this particular issue does not have the highest priority to me right now, but it’s good to know that the issue is recognized as a bug.


(Dženan Zukić) #29

I am not eager to do it, I only wanted to make sure somebody puts it on their to-do list. Because the discussion’s conclusion was that it should be done, without a hint as who should do it :smiley:


(Niels Dekker) #30

OK, I’ll have a look at how to fix the issue by moving code from GetCircles() to AfterThreadedGenerateData, when I have time. But first I would like you to have a look at this new proposed patch of mine: “ENH: Improved HoughTransform2DCirclesImageFilter accuracy, using Math::Round”, http://review.source.kitware.com/#/c/22784


(Niels Dekker) #31

The documentation says that the filter finds circles, as in the example data:

However, the filter can also be used to detect discs (a disc being a region bounded by a circle), assuming that the intensity within the disc is higher than outside the disc. Actually I think it will find a high-intensity-disc more easily than a circle that only has high intensity pixels on the boundary, as the algorithm uses gradient vectors. For a high-intensity disc, the gradient vectors all point away from the disc center. For a high-intensity-circle-boundary, there are gradient vectors pointing away from the circle center, as well as gradient vectors toward the circle center.

Should there be a note added to the documentation that the filter finds discs as well?

By the way, I did already adapt the test code, to directly pass an image of “discs” (generated by ‘CreateCircle’) to the Hough filter, instead of passing images of circle boundaries. I did so by removing the gradient filtering + thresholding preprocessing steps from the test: http://review.source.kitware.com/#/c/22784/3/Modules/Filtering/ImageFeature/test/itkHoughTransform2DCirclesImageTest.cxx That’s OK, right? It appears to perform just fine (or better) without those preprocessing steps.


(Dženan Zukić) #32

Expansion and clarification of the documentation is always welcome! The process is same as for code modification (submit to Gerrit etc).


(Niels Dekker) #33

OK, here is my proposed expansion of the documentation: “DOC: Mentioned that HoughTransform2DCirclesImageFilter also finds discs”, http://review.source.kitware.com/#/c/22821/

By the way, the link from the documentation main page of HoughTransform2DCirclesImageFilter to its example documentation page is broken. At https://itk.org/Doxygen/html/classitk_1_1HoughTransform2DCirclesImageFilter.html (at the section “Wiki Examples”).

If you can fix it, please do :slight_smile:


(Dženan Zukić) #34

@matt.mccormick has been updating documentation system recently. Matt, can you fix the broken link?


(Matt McCormick) #35

Here’s a patch to fix the broken documentation link.


(Niels Dekker) #36

FYI, I just had a closer look at how the patch I proposed last week, “PERF: Improved speed of copying and resizing NeighborhoodAllocator” http://review.source.kitware.com/#/c/22862/ would affect the performance of ITK’s Hough Transform.

I observed a 10% reduction (approximately) of the duration of filter->Update(), typically going from 50 to 45 seconds, on the following test program, using VS2017, Release, 64-bit:

#include <itkHoughTransform2DCirclesImageFilter.h>
#include <iostream>
#include <ctime>

int main()
{
  typedef unsigned char PixelType;
  const auto image = itk::Image<PixelType>::New();
  enum { sizeX = 4096, sizeY = sizeX };
  image->SetRegions({ sizeX, sizeY });
  image->Allocate();
  image->FillBuffer(1);

  const auto filter =
    itk::HoughTransform2DCirclesImageFilter<PixelType, unsigned char, double>::New();
  filter->SetInput(image);

  const auto time1 = std::time(nullptr);
  filter->Update();
  const auto time2 = std::time(nullptr);

  std::cout << "Duration: " << (time2 - time1) << " seconds" << std::endl;
}

HoughTransform2DCirclesImageFilter::GenerateData() indirectly does a lot of Neighborhood allocations, by calling GaussianDerivativeImageFunction::EvaluateAtIndex on each input pixel > Threshold.


Removal of `virtual` keywords from ConstNeighborhoodIterator
(Dženan Zukić) #37

As that patch does not make big changes, I am in favor of merging it before release. @matt.mccormick, @blowekamp, @fbudin?


(Matt McCormick) #38

This is a good patch, but we should keep to the feature freeze. Small patches still have unintended consequences. We want a robust release. It will be merged soon, and the next release will come soon enough.


(Niels Dekker) #39

GaussianDerivativeImageFunction::EvaluateAtIndex is called many times during an Update() of HoughTransform2DCirclesImageFilter. Can someone please explain this function to me?

Specifically, I see the center of all of its kernels is set to zero (both on line 267 and 277):

https://github.com/Kitware/ITK/blob/v4.13rc02/Modules/Core/ImageFunction/include/itkGaussianDerivativeImageFunction.hxx#L267

Do you know why it does so?

Update (19 December 2017): Do you agree that ‘centerval’ is always zero, after the first EvaluateAtIndex call? If so, doesn’t the following patch basically produce the same result as the original code (while being much faster)? “WIP: Removed blurring from GaussianDerivativeImageFunction::EvaluateAtIndex”, http://review.source.kitware.com/#/c/22926/


(Niels Dekker) #40

@dzenanz Thanks for recognizing the GaussianDerivativeImageFunction::EvaluateAtIndex bug, at http://review.source.kitware.com/#/c/22926/ Hope someone can still explain the intention of the function, and how to properly fix it!

There’s another major bug in GaussianDerivativeImageFunction, which directly affects HoughTransform2DCirclesImageFilter. Its protected helper function RecomputeGaussianKernel() may or may not use the image spacing. You see, it has a member m_UseImageSpacing, which is always true. It cannot be set by the user! That would suggest that the spacing is always used, right? RecomputeGaussianKernel() is called in the default-constructor of GaussianDerivativeImageFunction. But within the default-constructor, the image spacing is not yet known! So then the image spacing is still ignored! GaussianDerivativeImageFunction::SetSigma does call RecomputeGaussianKernel() again, but only when the new sigma is different from the previous one. When HoughTransform2DCirclesImageFilter does set sigma to 1.0 (which is the default), SetSigma usually does not call RecomputeGaussianKernel(), because sigma has not changed.

Note that SetSigma should be called after SetInputImage, because otherwise, RecomputeGaussianKernel() will use the spacing of the previous input image (if any). Note also that SetInputImage does not call RecomputeGaussianKernel().

See https://github.com/Kitware/ITK/blob/v4.13.0/Modules/Core/ImageFunction/include/itkGaussianDerivativeImageFunction.hxx