Interpolating values on a 2D plot when NOT a function (i.e. a looped shape)

Hi maestros,

I have a shape in x,y, that loops back on itself and meanders (without intersecting itself). An example would be a circle, a heart shape or the top view of a race car track.

I wish to find the values of y (usually there would be more than 1) for a given x value (or set of x values). I am aware that the spline function could help me if this was a function. I am also aware that a parameter 't' could be made to find (x,y,) in terms of 't', but don't see how this could help me.

Here is an example set of points in x,y if you want it:

x= [-0.2073 0.5721 1.4722 2.2293 2.6327 3.0899 3.7784 4.5382 5.5277 6.9211 8.3167 9.2896 9.0604 7.8213 6.6720 5.9258 5.1584 4.3335 3.8667 4.1657 5.3247 5.9008 5.0912 3.0462 0.6883 -0.9153 -1.9870 -3.2646 -4.0178 -4.6815 -6.2504 -7.6611 -8.3864 -8.6373 -8.3692 -7.5111 -6.5901 -6.3305 -6.6416 -7.2034 -7.1405 -6.1491 -4.8732 -3.8070 -3.0348 -2.4280 -1.5272 -0.8534 -0.4118 0.2505]

y = [ -2.3595 -2.9712 -3.0852 -2.7334 -2.2111 -1.6185 -1.3093 -1.3964 -1.6764 -1.9983 -1.8301 -0.9510 0.0921 0.7508 0.8846 0.5857 0.2114 0.1446 0.4953 1.1859 1.8096 2.5017 2.8755 2.7823 2.3773 1.7967 1.2743 1.2158 1.7293 2.3063 2.5951 2.4283 1.8694 1.1793 0.6032 0.2349 -0.0044 -0.3900 -1.1001 -1.7027 -2.3185 -2.7130 -2.6429 -2.1655 -1.6388 -1.4356 -1.4135 -1.6123 -2.1106 -2.7784]

Any help would be greatly appreciated, as I have spent a LONG time searching.

Mike

Answer 1

[Image Analyst]

Mike: Do you want something like this, where I parameterized x and y and fit a spline to both?

% function SmoothSplineCurve()
clc;  % Clear command window.
clear;  % Delete all variables.
close all;  % Close all figure windows except those 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.

[Michael]

Hi Image Analyst, thanks for your response but unfortunately I can't see it helping. I have actually made a very similar function to yours in the past! I am working with curves created from such a UI tool.

The resultant x-y coords before and after your "splineXY = spline..." line are essentially random as they fall wherever your chosen 'finerSpacing' paramater places them. The problem is that I want specific x values, not specific parameter values.

Following your code with something like:

yNew = spline(x,y,desiredXcoordinates); %(I have made up these variable names)

will return an error as there is not one-to-one mapping between x and y for looped shapes.

[Image Analyst]

I'm not sure I follow you. Why can't you use random points in spline()??? Just because I evenly spaced my output "x" values (which I called finerSpacing) doesn't mean that you have to. You can pass in whatever coordinates you want for that. They can be randomly spaced as far as I can tell. Did you even try it?

[Shubham Burde]

Hello Image Analyst, Which kind of spline is formed using the spline command which you used in the above program? splineXY = spline(originalSpacing, knots, finerSpacing); i mean how can i know whether it is quadratic, cubic or any other order?

[Image Analyst]

By looking in the help. It says "Cubic spline data interpolation" as the description, so that should answer your question.

[Ryan Spencer]

Hello Image Analyst. Fantastic code! Could you describe how to superimpose an image as the background?You can then use your code to trace/outline a geometry within the image. Thank you!

[Image Analyst]

Not sure what you mean. Start a new question and give an example of your inputs and your desired output.

Answer 2

[Bjorn Gustavsson]

Pseudo-coded stab at this:

   1, find all indices for which x(i1) < x_current and x(i1+1) >  x_current
      and x(i1) > x_current and x(i1+1) <  x_current and x(i1) == x_current
   2, Interpolate between the point-pairs

HTH

Answer 3

[gaoliming gao]

You may want to refer to this exchange file https://www.mathworks.com/matlabcentral/fileexchange/34874-interparc