B-spline fitting to 2D discrete data points (pixels of contour image)

I am trying to fit a B-spline to a set of ordered discrete data points which represent pixels of a contour extracted from a image.

The below code works fine for some simple shapes, but not for others (please see image here http://imgur.com/5IOVmEt

for examples). Why does this happen, and what would be a better way to approach this problem?

I am quite new to differential geometry and MATLAB in general, appreciate any insights or inputs. Thanks.

x = data(:, 1); y = data(:, 2); % data contains two columns representing x,y coordinates
plot(x, y, 'bo');
fittedmodel = fit(x, y, 'cubicinterp');
plot(fittedmodel, 'r-');

Answer 1

[Image Analyst]

Looks to me like the problem is you sorted your x coordinates and didn't keep them in the order that they were are you traverse around the contour.

Why are you using splines? Do you want to interpolate at a higher resolution and keep the curve smooth and make sure it still goes through each original point? If so, use splines. If you want to smooth the boundary so there are less sharp jaggy edges and spikes and don't care if it goes exactly through each original point, then you can do other things, like perhaps a Savitzky-Golay filter (the sgolay function in the Signal Processing Toolbox).

[Pradeep Kumar]

Hi Image Analyst,

I want to use splines, as I want the curve to pass through each of the discrete points and recreate a resolution independent version of the discrete contour. I do this so that I can later calculate the curvature and use it as a feature for matching.

The data points are sorted correctly to traverse the contour from one endpoint to the to other, I just double-checked.

Thanks!

[Image Analyst]

Your curve is not parameterized or else you would not be bouncing from the top curve to the bottom curve like you are. I have two spline demos that may help you.

Here's one spline demo that you can adapt:

% Demo to show spline interpolation.
% Clean up / initialize
clc; 
close all; 
clear all; 
workspace; % Display workspace panel.

% Create the original knot points.
x = 1:10;
y = rand (10,1);
% Plot it and show how the line has sharp bends.
plot(x, y, '-sr', 'LineWidth', 2);
set(gcf, 'Position', get(0,'Screensize')); % Maximize figure.

% Use splines to interpolate a smoother curve,
% with 10 times as many points,
% that goes exactly through the same data points.
samplingRateIncrease = 10;
newXSamplePoints = linspace(1, 10, 10 * samplingRateIncrease);
smoothedY = spline(x, y, newXSamplePoints);

% Plot smoothedY and show how the line is
% smooth, and has no sharp bends.
hold on; % Don't destroy the first curve we plotted.
plot(newXSamplePoints, smoothedY, '-ob');
title('Spline Interpolation Demo', 'FontSize', 20);

legend('Original Points', 'Spline Points');

% Mathworks Demo code from their Help
% x = 0:10;
% y = sin(x);
% xx = 0:.25:10;
% yy = spline(x,y,xx);
% plot(x,y,'o',xx,yy)

Here is the second one. You may like this one better:

% function SmoothSplineCurve()
clc;  % Clear command window.
clear;  % Delete all variables.
close all;  % Close all figure windows except those created by imtool.
imtool close all;  % Close all figure windows created by imtool.
workspace;  % Make sure the workspace panel is showing.
fontSize = 20;

% Change the current folder to the folder of this m-file.
if(~isdeployed)
        cd(fileparts(which(mfilename)));
end
set(gcf, 'Position', get(0,'Screensize')); % Enlarge figure to full screen.
set(gcf,'name','Spline Image Analysis Demo','numbertitle','off') 

axis([0 10 0 10])
xlabel('X', 'FontSize', fontSize);
ylabel('Y', 'FontSize', fontSize);
title('Spline Demo', 'FontSize', fontSize);
hold on
% Initially, the list of points is empty.
knots = [];
numberOfPointsClicked = 0;
% Prompt the user
message = sprintf('Left click to draw some vertex points.\nRight click the final point to finish drawing.');
uiwait(msgbox(message));
buttonThatWasClicked = 1;
% Enter a loop asking user to click on the knot vertexes.
while buttonThatWasClicked == 1
    [xKnot, yKnot, buttonThatWasClicked] = ginput(1);
    plot(xKnot, yKnot, 'ro', 'LineWidth', 2)
    numberOfPointsClicked = numberOfPointsClicked+1;
  % Make this coordinate a new column.
    knots(:, numberOfPointsClicked) = [xKnot; yKnot];
end

% Calculate the area within the blue spline curve.
% You do not need to connect the last point back to the first point.
x = knots(1, :);
y = knots(2, :);
areaOfPolygon = polyarea(x,y);

% Interpolate with a spline curve and finer spacing.
originalSpacing = 1 : numberOfPointsClicked;
% Make 9 points in between our original points that the user clicked on.
finerSpacing = 1 : 0.1 : numberOfPointsClicked;
% Do the spline interpolation.
splineXY = spline(originalSpacing, knots, finerSpacing);

% Plot the interpolated curve.
hold off;
plot(knots(1, :), knots(2, :), 'ro',...
  splineXY(1, :), splineXY(2, :), 'b+-', 'LineWidth', 2, 'MarkerSize', 16);
title('Blue Spline Between Red Knots', 'FontSize', fontSize);
legend('Knots', 'Spline');
xlabel('X', 'FontSize', fontSize);
ylabel('Y', 'FontSize', fontSize);
grid on;
hold off;
% Calculate the area within the blue spline curve.
% You do not need to connect the last point back to the first point.
x = splineXY(1, :);
y = splineXY(2, :);
areaInsideSplineCurve = polyarea(x,y);

% Give the area calculations.
message = sprintf('The area inside the polygon you drew is %.2f.\nThe area inside the blue spline curve is %.2f', ...
  areaOfPolygon, areaInsideSplineCurve);
fprintf(1, '%s', message); % Print to command window.
msgbox(message); % Show user via a popup message box.

[Pradeep Kumar]

A very comprehensive answer as always. The demos are cool, thanks! :)

[Evripidis]

This really helped me for a similar problem, as I am new to geometric computing. Is possible to have an example for the second demo using B-Splines and not cubic splines for interpolation/approximation?

Thanks a lot!

[Image Analyst]

I don't have a B spline demo. Never needed one before since the regular cubic spline always worked good enough.

[Evripidis]

Thanks for your quick reply! I suppose I will use cubic one because de Boor algorithm for b splines is too complex for me right now

Answer 2

[korakuti hanumanth]

Hi community,

I have two vectors and I would like to fit a B spline curve to:

y=[18.93000031 19.42000008 19.51000023 19.67000008 19.68000031 19.71999931];

x=[58.61111111 67.32055556 70.56194444 74.22694444 78.39388889 85.11555556];

I would also like to save the parameters of the fit at every point to use in regression analysis and also explain how to find the coefficients of b spline curve and some error techniques to satisfie the fitting.