I have a discrete contour represented by a set of points. The contour looks like a polygon but if you zoom you see that the edges are rugged (that's because it was obtained while working on a finite differences grid)
I would like to fit a polygon to this contour (which is closest in the least-squares sense). Are there any simple known algorithms to do this?
The Hough Transform is an image processing algorithm for extracting features for an image. The classical version of the algorithm is designed to extract lines from a binary image (such as this).
Given the ability to do this, you can make a script fairly easily. Here's one using the hough transform (and associated utilities) in Matlab. It probably requires the image processing toolbox, but hopefully you can track that down.
The only shortcut present in this code is that I told the "houghpeaks" function to search for 4 lines. You can be smarter than this by doing other peak detection methods in the Hough space (read the wiki if you're not familiar with the output of the Hough transform, https://en.wikipedia.org/wiki/Hough_transform).
One caveat is that the Hough transform will only give you whole lines. You'll have to post-process them to get line segments. This example works particularly well without any of that since the input was a convex polygon whose vertices were at/near the edge of the image.
clear all
close all
clc
% Load image (I suspect your input may be different, but the same ideas hold)
img = imread('polygon.png');
% Make a black/white image using the "blue" as the new white values
bw = img(:,:,1) ~= 255;
% Use hough transform to find lines
[H,theta,rho] = hough(bw);
P = houghpeaks(H,4);
lines = houghlines(bw,theta,rho,P);
% Visualization
imagesc(bw); colormap gray;
hold on;
for l = 1:length(lines)
plot([lines(l).point1(1) lines(l).point2(1)], ...
[lines(l).point1(2) lines(l).point2(2)],'r');
end
I am quite surprised why no-one mentioned the famous Douglas Peucker Algorithm for polyline simplification. Since you have contour points in hand, you could benefit from it directly. Contour approximation in OpenCV uses this method. See this for usage and you could also find a MATLAB implementation here or here.
There is a approach that uses notions from Discrete Geometry: Discrete Geometry is a discipline that works with objects defined as sets of pixels that try to mimic their standard counterparts. It defines discrete segments, discrete circles, discrete planes etc... In your case, there is an algorithm [1] that has a definition of what a discrete segment is, and that reconstructs in the input image the set of maximal segments, i.e. segments that cannot be further extended by adding new pixels to them. See also the extensions to fuzzy data [2,3]. The approach was successfully applied to 3D reconstruction from multiview images [4]
[1] I. DEBLED-RENNESSON, J.-P. REVEILLES, A linear algorithm for segmentation of digital curves, International Journal of Pattern Recognition and Artificial Intelligence, Volume 9, N. 6, December 1995.
[2] I. DEBLED-RENNESSON, J.-L. REMY and J. ROUYER-DEGLI, Linear segmentation of discrete curves into fuzzy segments, Discrete Applied Mathematics, 151:122-137, October 2005.
[3] I. DEBLED-RENNESSON, F. FESCHET and J. ROUYER-DEGLI, Optimal blurred segments decomposition of noisy shapes in linear time, Computers & Graphics, 30(1), 2006.
I wouldn't call it simple, but this paper offers a least squares fitting of data to polygons. If you're looking for accuracy, this will likely be the best option.
https://www.researchgate.net/publication/303329343_Least-squares_Fitting_of_Polygons
Active shape model is a slick algorithm for this kind of scenario.
Another approach is to break it up into 4 individual linear regression problems. You have sharp corners, so create an array of points, apply a basic smooth algorithm, and look at the slope of each point to its neighbor (ie derivative). Or if you didn't have sharp corners, fit a catmull-rom spline (boost c++ library has an implementation of this) and iterate over the spline looking at the derivative to determine the split points.
If you're looking for simple and fast, I agree with the Douglas Peuker Algorithm, but it needs to be tuned for the data set you are trying to fit, otherwise it might pick up corner points on the wrong line segments and hence allow unwanted error. You could use an iterative approach using multiple thresholds, measure the error, and try finer resolution iterations between points with the least error.
The hough transform and convex hull/convex defect have easy to use implementations in OpenCV and might also be helpful, as mentioned in the comment.
