ITK  6.0.0
Insight Toolkit
Examples/Segmentation/HoughTransform2DLinesImageFilter.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{HoughTransform2DLinesImageFilter} to find straight lines in a
// 2-dimensional image.
//
// First, we include the header files of the filter.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
// Software Guide : EndCodeSnippet
int
main(int argc, char * argv[])
{
if (argc < 4)
{
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 << " numberOfLines " << 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 = 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;
}
ImageType::Pointer localImage = reader->GetOutput();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// Once the image is loaded, we apply a
// \doxygen{GradientMagnitudeImageFilter} to segment edges. This casts
// the input image using a \doxygen{CastImageFilter}.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
using CastingFilterType =
auto caster = CastingFilterType::New();
std::cout << "Applying gradient magnitude filter" << std::endl;
using GradientFilterType =
itk::GradientMagnitudeImageFilter<AccumulatorImageType,
AccumulatorImageType>;
auto gradFilter = GradientFilterType::New();
caster->SetInput(localImage);
gradFilter->SetInput(caster->GetOutput());
gradFilter->Update();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// The next step is to apply a threshold filter on the gradient magnitude
// image to keep only bright values. Only pixels with a high value will be
// used by the Hough transform filter.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
std::cout << "Thresholding" << std::endl;
auto threshFilter = ThresholdFilterType::New();
threshFilter->SetInput(gradFilter->GetOutput());
threshFilter->SetOutsideValue(0);
unsigned char threshBelow = 0;
unsigned char threshAbove = 255;
threshFilter->ThresholdOutside(threshBelow, threshAbove);
threshFilter->Update();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We create the HoughTransform2DLinesImageFilter 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>;
auto houghFilter = HoughTransformFilterType::New();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We set the input to the filter to be the output of the
// ThresholdImageFilter. We set also the number of lines 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. SetDiscRadius() defines the radius of this disc.
//
// The output of the filter is the accumulator.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
houghFilter->SetInput(threshFilter->GetOutput());
houghFilter->SetNumberOfLines(std::stoi(argv[3]));
if (argc > 4)
{
houghFilter->SetVariance(std::stod(argv[4]));
}
if (argc > 5)
{
houghFilter->SetDiscRadius(std::stod(argv[5]));
}
houghFilter->Update();
AccumulatorImageType::Pointer localAccumulator = houghFilter->GetOutput();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We can also get the lines as \doxygen{LineSpatialObject}. The
// \code{GetLines()} function return a list of those.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
HoughTransformFilterType::LinesListType lines;
lines = houghFilter->GetLines();
std::cout << "Found " << lines.size() << " line(s)." << std::endl;
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We can then allocate an image to draw the resulting lines as binary
// objects.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
using OutputPixelType = unsigned char;
using OutputImageType = itk::Image<OutputPixelType, Dimension>;
auto localOutputImage = OutputImageType::New();
OutputImageType::RegionType region(localImage->GetLargestPossibleRegion());
localOutputImage->SetRegions(region);
localOutputImage->CopyInformation(localImage);
localOutputImage->Allocate(true); // initialize buffer to zero
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We iterate through the list of lines and we draw them.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
using LineIterator =
HoughTransformFilterType::LinesListType::const_iterator;
LineIterator itLines = lines.begin();
while (itLines != lines.end())
{
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We get the list of points which consists of two points to represent a
// straight line. Then, from these two points, we compute a fixed point
// $u$ and a unit vector $\vec{v}$ to parameterize the line.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
using LinePointListType =
HoughTransformFilterType::LineType::LinePointListType;
LinePointListType pointsList = (*itLines)->GetPoints();
LinePointListType::const_iterator itPoints = pointsList.begin();
double u[2];
u[0] = (*itPoints).GetPositionInObjectSpace()[0];
u[1] = (*itPoints).GetPositionInObjectSpace()[1];
itPoints++;
double v[2];
v[0] = u[0] - (*itPoints).GetPositionInObjectSpace()[0];
v[1] = u[1] - (*itPoints).GetPositionInObjectSpace()[1];
double norm = std::sqrt(v[0] * v[0] + v[1] * v[1]);
v[0] /= norm;
v[1] /= norm;
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We draw a white pixels in the output image to represent the line.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
itk::Size<2> size =
localOutputImage->GetLargestPossibleRegion().GetSize();
float diag =
std::sqrt(static_cast<float>(size[0] * size[0] + size[1] * size[1]));
for (auto i = static_cast<int>(-diag); i < static_cast<int>(diag); ++i)
{
localIndex[0] = static_cast<long>(u[0] + i * v[0]);
localIndex[1] = static_cast<long>(u[1] + i * v[1]);
localOutputImage->GetLargestPossibleRegion();
if (outputRegion.IsInside(localIndex))
{
localOutputImage->SetPixel(localIndex, 255);
}
}
itLines++;
}
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We setup a writer to write out the binary image created.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
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::CastImageFilter
Casts input pixels to output pixel type.
Definition: itkCastImageFilter.h:100
itk::GradientMagnitudeImageFilter
Computes the gradient magnitude of an image region at each pixel.
Definition: itkGradientMagnitudeImageFilter.h:42
itk::Size< 2 >
itkImageFileReader.h
itkMinimumMaximumImageCalculator.h
itkImageRegionIterator.h
itkCastImageFilter.h
itkThresholdImageFilter.h
itk::ImageFileReader
Data source that reads image data from a single file.
Definition: itkImageFileReader.h:75
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
itk::HoughTransform2DLinesImageFilter
Performs the Hough Transform to find 2D straight lines in a 2D image.
Definition: itkHoughTransform2DLinesImageFilter.h:62
itk::ThresholdImageFilter
Set image values to a user-specified value if they are below, above, or outside threshold values.
Definition: itkThresholdImageFilter.h:71
itkImageFileWriter.h
itkHoughTransform2DLinesImageFilter.h
itk::Image
Templated n-dimensional image class.
Definition: itkImage.h:88
New
static Pointer New()
itk::GTest::TypedefsAndConstructors::Dimension2::Dimension
constexpr unsigned int Dimension
Definition: itkGTestTypedefsAndConstructors.h:44
itkDiscreteGaussianImageFilter.h
itk::Size::GetSize
const SizeValueType * GetSize() const
Definition: itkSize.h:169
itkGradientMagnitudeImageFilter.h