# 3D ellipse hough transform

**URL:** <https://discourse.itk.org/t/3d-ellipse-hough-transform/2861>\
**Category:** Algorithms\
**Created:** [March 20, 2020, 6:38pm UTC](https://discourse.itk.org/t/3d-ellipse-hough-transform/2861 "2020-03-20T18:38:15Z")\
**Posts on this page:** 2\
**Page:** 1

<div class="post-metadata">

**Author:** ![basal](https://discourse.itk.org/letter_avatar_proxy/v4/letter/b/a8b319/32.png) [@basal](https://discourse.itk.org/u/basal)\
**Post date:** [March 20, 2020, 6:38pm UTC](https://discourse.itk.org/t/3d-ellipse-hough-transform/2861/1 "2020-03-20T18:38:16Z")

</div>

Hi all, I want to extend 2D hough transform to 3D and also use ellipse function instead of circles in ITK based on c++ to detect ellipse shapes in my image. I changed some lines of the 2D Hough transform code for ellipse detection but it is not working and I am getting errors. Can you help me with that?  
Moreover, I do not have any idea about how to change 2D Hough transform to 3D. I really need these two changes and my time is so limited. I appreciate it if you could help me with these issues.  
Below is the 2D hough transform code I edited to detect ellipse but it’s not working.

```auto
// 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
#include "itkHoughTransform2DCirclesImageFilter.h"
// Software Guide : EndCodeSnippet
#include "itkImageFileReader.h"
#include "itkEllipseSpatialObject.h"
#include "itkImageFileWriter.h"
#include "itkImageRegionIterator.h"
#include "itkThresholdImageFilter.h"
#include "itkMinimumMaximumImageCalculator.h"
#include "itkGradientMagnitudeImageFilter.h"
#include "itkDiscreteGaussianImageFilter.h"
#include <list>
#include "itkCastImageFilter.h"
#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;
	using RadiusPixelType = float;
	constexpr unsigned int Dimension = 2;
	using ImageType = itk::Image< PixelType, Dimension >;
	ImageType::IndexType localIndex;
	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 >;
	ReaderType::Pointer reader = ReaderType::New();
	reader->SetFileName(argv[1]);
	try
	{
		reader->Update();
	}
	catch (itk::ExceptionObject & excep)
	{
		std::cerr << "Exception caught !" << std::endl;
		std::cerr << excep << std::endl;
		return EXIT_FAILURE;
	}
	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 =
		itk::HoughTransform2DCirclesImageFilter<PixelType,
		AccumulatorPixelType,
		RadiusPixelType>;
	HoughTransformFilterType::Pointer 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();
	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 EllipseSpatialObject<2>;
	EllipseSpatialObject<2> = houghFilter->GetCircles();
	std::cout << "Found " << EllipseSpatialObject<2>.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 >;
	OutputImageType::Pointer localOutputImage = OutputImageType::New();
	OutputImageType::RegionType region;
	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 != EllipseSpatialObject<2>.end())
	{
		std::cout << "Center: ";
		std::cout << (*itCircles)->GetObjectToParentTransform()->GetOffset()
			<< std::endl;
		std::cout << "Radius: " << (*itCircles)->GetRadius()[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 <= 2 * vnl_math::pi; angle += vnl_math::pi / 60.0)
		{
			localIndex[0] =
				(long int)((*itCircles)->GetObjectToParentTransform()->GetOffset()[0]
					+ (*itCircles)->GetRadius()[0] * cos(angle));
			localIndex[1] =
				(long int)((*itCircles)->GetObjectToParentTransform()->GetOffset()[1]
					+ (*itCircles)->GetRadius()[0] * sin(angle));
			OutputImageType::RegionType region =
				localOutputImage->GetLargestPossibleRegion();
			if (region.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 >;
	WriterType::Pointer writer = WriterType::New();
	writer->SetFileName(argv[2]);
	writer->SetInput(localOutputImage);
	try
	{
		writer->Update();
	}
	catch (itk::ExceptionObject & excep)
	{
		std::cerr << "Exception caught !" << std::endl;
		std::cerr << excep << std::endl;
		return EXIT_FAILURE;
	}
	// Software Guide : EndCodeSnippet
	return EXIT_SUCCESS;
}

```

---

<div class="post-metadata">

**Author:** ![beatrice](https://discourse.itk.org/user_avatar/discourse.itk.org/beatrice/32/2991_2.png) [@beatrice](https://discourse.itk.org/u/beatrice)\
**Post date:** [October 12, 2022, 10:36am UTC](https://discourse.itk.org/t/3d-ellipse-hough-transform/2861/2 "2022-10-12T10:36:06Z")

</div>

Hello, did you finally manage to solve the problem? Thanks
