Fit Curve to a "C" shape

Hey all,

I'm trying to use the curve fitting utilities to fit data although I'm having trouble. To make a simple example I am trying to fit a "C" shape with the code below. As you can see in the attached figure, the fitted curve does not follow the points, I'm assuming because such a curve would return two values for f(x) on the interval of x equals -1 to 1. My intended workaround is to make two splines, but this quickly becomes annoying for my actual problem, so I'm looking for a solution where Matlab draws a curve through the points in the order that are entered (kind of like the plotting utility). Any thoughts would be appreciated.

Thanks,

Chris

x = [1 0 -1 0 1];
y = [-1 -1 0 1 1];

f = fit(x',y','smoothingspline');

plot(f,x,y,'o')
xlim([-2 2])
ylim([-2 2])

[Christopher Axten]

Found an answer online for anyone else with the issue:

https://www.mathworks.com/matlabcentral/answers/411231-connecting-data-points-in-a-smooth-curve

They really should just add an example like this to the spline documentation page for clarity.

Answer 1

[John D'Errico]

A spline produces a FUNCTION, SINGLE VALUED. That is, for any x, it is assumed there is a single y. Is that the case for you? NO.

In that case, you need to convert this to a problem that is so. Or you need to use a tool designed to handle non-single valued relationships, so general curvilinear arcs in the plane.

x = [1 0 -1 0 1];
y = [-1 -1 0 1 1];

[theta,R] = cart2pol(x,y);
theta = unwrap(theta);

We need to be careful, and use unwrap here. It will resolve jumps in polar angle as the curve wraps around. Otherwise, spline will get confused.

The simple solution, since your function already lives nicely in a non-translated polar coordinate frame around the origin, is to do as I did above, now just use an interpolant. Spline will suffice, or we could be even more lazy and use interp1.

r_of_theta = spline(theta,R);
th_interp = linspace(min(theta),max(theta),100);
[x_interp, y_interp] = pol2cart(th_interp,fnval(r_of_theta,th_interp));
plot(x,y,'ro',x_interp,y_interp,'b-')
axis equal
grid on

Remember to use unwrap!

This works nicely enough, at least for this simple problem. For more complicated problems, you may want to try my interparc code, as found on the file exchange. It could handle more complicated problems, even ones where the curve crosses itself. Find it here:

https://www.mathworks.com/matlabcentral/fileexchange/34874-interparc

xyint = interparc(100,x,y);
plot(x,y,'ro',xyint(:,1),xyint(:,2),'b-')
axis equal
grid on

It produced a subtly different curve shape, because the spline it uses is a subtly different spline model. But the nice thing about interparc is it does not ask you to do anything special like understand a conversion to polar coordinates, or to understand why unwrap was necessary in my first fit.

[Christopher Axten]

John,

When I was writing up the question I was going to include the part about a function requiring one value of f(x) for every x (remembering way back to pre-calc), so thank you for reaffiriming that's why Matlab is having an issue here.

The points I'm trying to fit actually look like this, and it is hard to tell but the points in the lower left (e.g. 0,.3; 0,0.6, etc.) are just to the right of a couple points in the upper left, which makes the fit spline zigzag between them. I found another approach, basically what Sulaymon outlined, but I like the converting to polar idea because it might handle the sharp change better (we'll see). I don't need to capture that spike perfectly, but the zigzagging I had before was not adequate, so I think your ideas will help with that.

Thank you!

Chris

[Walter Roberson]

Is that intended to be a near-vertical line that slopes slightly leftwards, that meets a curve that heads right?

Is it intended to be a continuous curve that has a knot?

Is it intended to be a continuous curve that goes something like

**
*  *
*    *
 *    *
 *     **
 *       ****
 *

that is, something that has a cap at the top left that meets both segments?

[Christopher Axten]

Hi Walter,

It is a continuous curve similar to the image you made. I was able to essentially what I want with the same approach Sulaymon mentioned. A zoomed in version with the spline fitted is below. I don't like the flat part at the peak, but I probably just need to add more points to round it. I am thinking using polar coordinates with linspace would add points where the curvature is high, so that is what I am about to try. I have not seen a tool to, for example, extract 100 points on the spline but bias them to capture the curvature as best as possible, but I may try to look for something like that as well.

[Christopher Axten]

Using polar coordinates is not working as well as I hoped, but I think the spline I get in the above image will do the job.

[John D'Errico]

As always, the real problem behind the question is wildly different from the question asked.

What you have not done is provide the data. But it appears you still have a multi-valued function, so spline cannot fit it.

Also, NO. You do not want to use a polar coordinate transformation. I would try interparc, or at least an interpolant based on the methods used by interparc.

If you post the data, I can see what method may work reasonably well.

[Image Analyst]

Like you John, I asked for his actual data also. He answered me without posting it so evidently he does not want to post it, so I guess we're all done here and he'll continue on his own.

[Christopher Axten]

Hi John and Image?,

I have attached my current code and data here if you are curious. I was slow to share because I pretty quickly got the spline I have above and it addressed my need in regard to the original question. I am now trying to differentiate that spline, but I didn't include this in the question so I didn't bring it up. I guess I should have included my whole problem here, but I figured I would get quicker, simplier responses with a clear and simple example like the "C" above.

I tried using polar coord, including shifting the coordinates, which was a couple minutes wasted but hey I learned something.

The curve and its derivative are what I ultimately need, which seemed simple in Matlab if have a spline but doesn't seem to working as intended so I am in the process of figuring it out.

[Image Analyst]

This doesn't look like a C shape. Anyway, are you trying to do

1. a fit (regression where the points may not necessarily go through your training points), or
2. an interpolation (where you have additional points in between training points but also goes through the training points as well)?

Which one? I'm betting on interpolation. You have a flat-topped red curve because you're interpolating at points that do not include your training points. I don't have the curve fitting toolbox so I can't run your code. What x values do you want to interpolate? Some set of x values that does not include the training points, or additional x values in addition to the training points? Seems like the latter since you didn't like the flat red top. So I think you can do

xDesired = sort(unique([xTraining, xAdditional]));
yDesired = spline(x, y, xDesired);

[Christopher Axten]

You're right it doesn't, my intention with the "C" was to find a solution for splining something that isn't a fuction, like in the lower x coordinates of the data, so in the future I'll be more forward with the exact problem at hand.

I am ultimately looking for the derivative of the data, so I was looking for a fit that I could differentiate. I initially used a second order finite difference to approximate it (pasted below) but it doesn't do very well near the front, as you would imagine. I don't know how much better a spline can do, but at least it would be smoother. There will be a point where the derivative (dUe/dx) is infinite because of how x doubles back, but I can avoid that point.

NOTE: I did figure out final answer to what I'm trying to do. It is unlikely your background is in the aerospace industry, but I was trying to work with the edge velocity from a stagnation point on an airfoil. The work around is to start x at the stagnation point, rather than the geometric leading edge, which makes the whole thing a function and all of the difficult parts go away.

Thank you to John, Image, and Walter for your help and thoughts!

dUdx = zeros(1,length(x));
for d = 2:length(x)-1
    dUdx(d) = (Ue(d+1)-Ue(d-1))/(2*x(d+1)-x(d-1));
end

Answer 2

[Walter Roberson]

You will not be able to do this in cartesian coordinates. You will need to transform to some other coordinate system, such as angle (independent) and magnitude (dependent) relative to the centroid.

Answer 3

[Sulaymon Eshkabilov]

Just a piecewise interpolation can be done in this case, but not fit model calculation:

clearvars
x = [1 0 -1 0 1];
y = [-1 -1 0 1 1];
N = 1:numel(x);
xy = [x;y];
SP = makima(N, xy);
%SP = pchip(N, xy);
IN = linspace(1,numel(x));
INTER = ppval(SP, IN);
plot(x,y,'bd', INTER(1,:),INTER(2,:), 'r-', 'linewidth', 2)
legend('DATA','Interpolation', 'location', 'best')

[Christopher Axten]

Sulaymon,

Thank you for the response! This is pretty close to what I found in another question after I posted my question, and it worked pretty well. My actual problem is more complex than the example "C" I presented, so I am going to compare this approach and the polar transformation approach that John outlined.

Thanks,

Chris

Answer 4

[Image Analyst]

You can fit the data to an ellipse, like in my attached demo. Then when you construct the curve you can just crop off any fitted values with x > your max x. Would an ellipse work for you?

Please attach a .mat file with your actual (x,y) data - sometimes simple data like you posted are too "perfect" and cause problems, like if your "simple" x or y values are perfect duplicates but your actual data don't have duplicates.

[Christopher Axten]

Hey, thank you for the response. I posted an image of the points I am trying to fit under John's comment. I don't think an ellipse would fit my actual problem well, but it would probably work great for the example I based the question on. Sorry for the miscommunication there, I was trying to keep it simple for brevity.