# Curve to approximate a sorted set of points

• Jun 1st 2013, 01:39 PM
mimoio
Curve to approximate a sorted set of points
Hi everyone,

i have a sorted list of N 2D points: {(x0,y0), (x1,y1), ... (xN,yN)}

and i would like to find a parametrized curve c(t) = (x(t), y(t)) that approximates the points following the preset order.

The curve may not be a function, that is, it may have y^{-1} (y_i) = y^{-1} (y_j) = x_k, points may share the same x values.

I have been looking for splines and other curves but i don't really what solutions have already been implemented and are robust and well known.

Is there something similar with 3D points?

• Jun 2nd 2013, 01:37 AM
mimoio
Re: Curve to approximate a sorted set of points
Sorry i didn't mean to approximate the points, but to PASS through them. Thanks and apologies.
• Jun 3rd 2013, 07:00 AM
johng
Re: Curve to approximate a sorted set of points
Hi,
I'd say the most common solution to your interpolation problem is to use splines. The simplest case is that of "natural" cubic splines. Eric Weisstein's page Cubic Spline -- from Wolfram MathWorld has all the necessary information. The solution there shows interpolation of y coordinates only. You'll want to also produce the spline for the x coordinates.

I've attached a drawing that uses 5 points -- the curve is closed in the sense that a 6th point, equal to the first point, is used in the interpolation.

Attachment 28516
• Jun 4th 2013, 10:59 AM
mimoio
Re: Curve to approximate a sorted set of points
Hi ahd thx jonhg,

this seems to be what i am looking for. I'll try to work on the algorithm in your wolfram link.

I've also had a look at matlab function "spline". Do you think that it fits my purpose? I ask this cause i made some tests and it seemed not to accout for the sorted points.
• Jun 4th 2013, 02:41 PM
mimoio
Re: Curve to approximate a sorted set of points
hi, finally i coded the algorithm u showed me... and it works fantastically :D!!! Thx a lot!

Attachment 28523

I paste the matlab code for anyone to use. Its divided in two files:

First file, where to execute it from, u can name it as u want:
----------------------------------------------------------------------------------
close all;
clear all;

plot_step = 0.1;

points = [
1,0;
1,1;
2,3;
4,5;
5,3;
4,3;
4,4;
3,4
];

plot(points(:,1),points(:,2),'or')
hold on;

px = points(:,1)';
py = points(:,2)';

[curve_params_x, curve_points_x] = compute_onedimension_cubic_spline(px, plot_step);
[curve_params_y, curve_points_y] = compute_onedimension_cubic_spline(py, plot_step);

plot(curve_points_x,curve_points_y,'-b')
axis([min(px)-1 max(px)+1 min(py) max(py)+1])
hold on;
------------------------------------------ END FILE ---------------------------------------------

Second file, name it "compute_onedimension_cubic_spline.m"
---------------------------------------------------------------------------------------
function [curve_params, curve_points] = compute_onedimension_cubic_spline(data, step)

% INPUT DATA
% data = row
% step = step to compute the curve_points

n = size(data,2) - 1;
if n < 1
error('A minimum of 2 points is needed to compute a curve')
end

% Construct tridiagonal matrix
A = zeros(n+1);
A(1,1) = 2; A(1,2) = 1;
for i = 2:n
A(i,i-1) = 1; A(i,i) = 4; A(i,i+1) = 1;
end
A(n+1,n) = 1; A(n+1,n+1) = 2;

% Construct points difference vector
pd = zeros(n+1,1);
pd(1) = 3*(data(2) - data(1));
for i = 2:n
pd(i) = 3*(data(i+1) - data(i-1));
end
pd(n+1) = 3*(data(n+1) - data(n));

% Compute vector D
D = A\pd;

%[A, D, pd]

% Compute parameters of cubic pieces of lines: [a(i), b(i), c(i), d(i)]
curve_params = zeros(n,4);
for i=1:n
a = data(i);
b = D(i);
c = 3*(data(i+1)-data(i)) - 2*D(i) - D(i+1);
d = 2*(data(i)-data(i+1)) + D(i) + D(i+1);
curve_params(i,:) = [a, b, c, d];
end
%curve_params

% Curve points
n_steps = 1/step;
curve_points = zeros(1,n*n_steps+1);
for i=1:n
a = curve_params(i,1);
b = curve_params(i,2);
c = curve_params(i,3);
d = curve_params(i,4);
for j = 1:n_steps
t = (j-1)*step;
curve_points((i-1)*n_steps + j) = a + b*t + c*t*t + d*t*t*t;
end
end
curve_points(n*n_steps+1) = a + b + c + d;
%curve_points
-------------------------------------------- END FILE -------------------------------------------