ITK  6.0.0
Insight Toolkit
Examples/Segmentation/HoughTransform2DCirclesImageFilter.cxx
/*=========================================================================
*
* Copyright NumFOCUS
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* https://www.apache.org/licenses/LICENSE-2.0.txt
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*
*=========================================================================*/
// Software Guide : BeginLatex
//
// This example illustrates the use of the
// \doxygen{HoughTransform2DCirclesImageFilter} to find circles in a
// 2-dimensional image.
//
// First, we include the header files of the filter.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
// Software Guide : EndCodeSnippet
#include <list>
#include "itkMath.h"
int
main(int argc, char * argv[])
{
if (argc < 6)
{
std::cerr << "Missing Parameters " << std::endl;
std::cerr << "Usage: " << argv[0] << std::endl;
std::cerr << " inputImage " << std::endl;
std::cerr << " outputImage" << std::endl;
std::cerr << " numberOfCircles " << std::endl;
std::cerr << " radius Min " << std::endl;
std::cerr << " radius Max " << std::endl;
std::cerr << " sweep Angle (default = 0)" << std::endl;
std::cerr << " SigmaGradient (default = 1) " << std::endl;
std::cerr << " variance of the accumulator blurring (default = 5) "
<< std::endl;
std::cerr
<< " radius of the disk to remove from the accumulator (default = 10) "
<< std::endl;
return EXIT_FAILURE;
}
// Software Guide : BeginLatex
//
// Next, we declare the pixel type and image dimension and specify the
// image type to be used as input. We also specify the image type of the
// accumulator used in the Hough transform filter.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
using PixelType = unsigned char;
using AccumulatorPixelType = unsigned int;
using RadiusPixelType = float;
constexpr unsigned int Dimension = 2;
using AccumulatorImageType = itk::Image<AccumulatorPixelType, Dimension>;
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We setup a reader to load the input image.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
using ReaderType = itk::ImageFileReader<ImageType>;
auto reader = ReaderType::New();
reader->SetFileName(argv[1]);
try
{
reader->Update();
}
catch (const itk::ExceptionObject & excep)
{
std::cerr << "Exception caught !" << std::endl;
std::cerr << excep << std::endl;
return EXIT_FAILURE;
}
const ImageType::Pointer localImage = reader->GetOutput();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We create the HoughTransform2DCirclesImageFilter based on the pixel
// type of the input image (the resulting image from the
// ThresholdImageFilter).
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
std::cout << "Computing Hough Map" << std::endl;
using HoughTransformFilterType =
AccumulatorPixelType,
RadiusPixelType>;
auto houghFilter = HoughTransformFilterType::New();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We set the input of the filter to be the output of the
// ImageFileReader. We set also the number of circles we are looking for.
// Basically, the filter computes the Hough map, blurs it using a certain
// variance and finds maxima in the Hough map. After a maximum is found,
// the local neighborhood, a circle, is removed from the Hough map.
// SetDiscRadiusRatio() defines the radius of this disc proportional to
// the radius of the disc found. The Hough map is computed by looking at
// the points above a certain threshold in the input image. Then, for each
// point, a Gaussian derivative function is computed to find the direction
// of the normal at that point. The standard deviation of the derivative
// function can be adjusted by SetSigmaGradient(). The accumulator is
// filled by drawing a line along the normal and the length of this line
// is defined by the minimum radius (SetMinimumRadius()) and the maximum
// radius (SetMaximumRadius()). Moreover, a sweep angle can be defined by
// SetSweepAngle() (default 0.0) to increase the accuracy of detection.
//
// The output of the filter is the accumulator.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
houghFilter->SetInput(reader->GetOutput());
houghFilter->SetNumberOfCircles(std::stoi(argv[3]));
houghFilter->SetMinimumRadius(std::stod(argv[4]));
houghFilter->SetMaximumRadius(std::stod(argv[5]));
if (argc > 6)
{
houghFilter->SetSweepAngle(std::stod(argv[6]));
}
if (argc > 7)
{
houghFilter->SetSigmaGradient(std::stoi(argv[7]));
}
if (argc > 8)
{
houghFilter->SetVariance(std::stod(argv[8]));
}
if (argc > 9)
{
houghFilter->SetDiscRadiusRatio(std::stod(argv[9]));
}
houghFilter->Update();
const AccumulatorImageType::Pointer localAccumulator =
houghFilter->GetOutput();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We can also get the circles as \doxygen{EllipseSpatialObject}. The
// \code{GetCircles()} function return a list of those.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
HoughTransformFilterType::CirclesListType circles;
circles = houghFilter->GetCircles();
std::cout << "Found " << circles.size() << " circle(s)." << std::endl;
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We can then allocate an image to draw the resulting circles as binary
// objects.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
using OutputPixelType = unsigned char;
using OutputImageType = itk::Image<OutputPixelType, Dimension>;
auto localOutputImage = OutputImageType::New();
region.SetSize(localImage->GetLargestPossibleRegion().GetSize());
region.SetIndex(localImage->GetLargestPossibleRegion().GetIndex());
localOutputImage->SetRegions(region);
localOutputImage->SetOrigin(localImage->GetOrigin());
localOutputImage->SetSpacing(localImage->GetSpacing());
localOutputImage->Allocate(true); // initializes buffer to zero
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We iterate through the list of circles and we draw them.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
using CirclesListType = HoughTransformFilterType::CirclesListType;
CirclesListType::const_iterator itCircles = circles.begin();
while (itCircles != circles.end())
{
std::cout << "Center: ";
std::cout << (*itCircles)->GetCenterInObjectSpace() << std::endl;
std::cout << "Radius: " << (*itCircles)->GetRadiusInObjectSpace()[0]
<< std::endl;
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We draw white pixels in the output image to represent each circle.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
for (double angle = 0; angle <= itk::Math::twopi;
angle += itk::Math::pi / 60.0)
{
(*itCircles)->GetCenterInObjectSpace();
localIndex[0] = itk::Math::Round<IndexValueType>(
centerPoint[0] +
(*itCircles)->GetRadiusInObjectSpace()[0] * std::cos(angle));
localIndex[1] = itk::Math::Round<IndexValueType>(
centerPoint[1] +
(*itCircles)->GetRadiusInObjectSpace()[0] * std::sin(angle));
const OutputImageType::RegionType outputRegion =
localOutputImage->GetLargestPossibleRegion();
if (outputRegion.IsInside(localIndex))
{
localOutputImage->SetPixel(localIndex, 255);
}
}
itCircles++;
}
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We setup a writer to write out the binary image created.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
using WriterType = itk::ImageFileWriter<ImageType>;
auto writer = WriterType::New();
writer->SetFileName(argv[2]);
writer->SetInput(localOutputImage);
try
{
writer->Update();
}
catch (const itk::ExceptionObject & excep)
{
std::cerr << "Exception caught !" << std::endl;
std::cerr << excep << std::endl;
return EXIT_FAILURE;
}
// Software Guide : EndCodeSnippet
return EXIT_SUCCESS;
}
Pointer
SmartPointer< Self > Pointer
Definition: itkAddImageFilter.h:93
itk::GTest::TypedefsAndConstructors::Dimension2::PointType
ImageBaseType::PointType PointType
Definition: itkGTestTypedefsAndConstructors.h:51
itkImageFileReader.h
itkMinimumMaximumImageCalculator.h
itkImageRegionIterator.h
itkCastImageFilter.h
itkThresholdImageFilter.h
itk::IndexValueType
long IndexValueType
Definition: itkIntTypes.h:93
itk::ImageFileReader
Data source that reads image data from a single file.
Definition: itkImageFileReader.h:75
itkHoughTransform2DCirclesImageFilter.h
itk::GTest::TypedefsAndConstructors::Dimension2::IndexType
ImageBaseType::IndexType IndexType
Definition: itkGTestTypedefsAndConstructors.h:50
itk::ImageFileWriter
Writes image data to a single file.
Definition: itkImageFileWriter.h:90
itk::GTest::TypedefsAndConstructors::Dimension2::RegionType
ImageBaseType::RegionType RegionType
Definition: itkGTestTypedefsAndConstructors.h:54
itkImageFileWriter.h
itk::ExceptionObject
Standard exception handling object.
Definition: itkExceptionObject.h:50
itk::Image
Templated n-dimensional image class.
Definition: itkImage.h:88
itk::HoughTransform2DCirclesImageFilter
Performs the Hough Transform to find circles in a 2D image.
Definition: itkHoughTransform2DCirclesImageFilter.h:61
New
static Pointer New()
itk::ImageRegion::SetSize
void SetSize(const SizeType &size)
Definition: itkImageRegion.h:202
itk::Math::pi
static constexpr double pi
Definition: itkMath.h:66
itk::GTest::TypedefsAndConstructors::Dimension2::Dimension
constexpr unsigned int Dimension
Definition: itkGTestTypedefsAndConstructors.h:44
itkDiscreteGaussianImageFilter.h
itk::Math::twopi
static constexpr double twopi
Definition: itkMath.h:68
itkMath.h
itkGradientMagnitudeImageFilter.h