How to fit an equation (catenary equation in 2D) to a 3D dataset (experimental data of a catenary curve in 3D) ?

Question

VIGNESH BALAJI 2023-7-20

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1998593-how-to-fit-an-equation-catenary-equation-in-2d-to-a-3d-dataset-experimental-data-of-a-catenary-cu

评论： Matt J 2023-7-28

I have an experimental data (3D data) of a catenary curve (basically a rope with 2 attached ends and hanging freely, data got from markers on them).

I want to fit this data with a catenary equation y = a*cosh(x/a).

a is parameterized by an equation stating the length of the physical cable, (length/2) - a* sinh(x/a) = 0

I tried using curve fitting tool and it is asking me an equation in terms of x,y and z. If I give such an equation then it is a surface.

Further, I don't know how to give an equation to parameterize parameter a in the curve fitting tool.

My ultimate requirement is, I want to fit a line (2D equation) in 3D space.

Can anyone please let me know the right tool for this ?

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Matt J 2023-7-20

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1998593-how-to-fit-an-equation-catenary-equation-in-2d-to-a-3d-dataset-experimental-data-of-a-catenary-cu#answer_1276333

编辑：Matt J 2023-7-20

在 MATLAB Online 中打开

In the Examples tab of this FEX page,

https://www.mathworks.com/matlabcentral/fileexchange/87584-object-oriented-tools-for-fitting-conics-and-quadrics

there is a section Fitting a 2D shape to 3D Points. In particular these lines of code coould easily be adapted to your problem,

pFit=planarFit(XYZ0);%Preliminary plane fit
  xy0=pFit.project2D(XYZ0); %Map measured 3D samples to 2D
  
eFit=ellipticalFit(xy0); %Perform ellipse fit in 2D
  XYZ=pFit.unproject3D(  cell2mat(eFit.sample(1:360))  ); %Post-sample the ellipse fit and map back to 3D

except that instead of the ellipse fit, you would fit with your own model.

6 个评论
显示 4更早的评论隐藏 4更早的评论

VIGNESH BALAJI 2023-7-27

planar_fit_structure.png

Thank you @Matt J it worked. I was trying it with plot3 initially instead of plot and internally the structure of pFit is like this, like the image attached below. How does it plot using this structure, can you explain me how.

Matt J 2023-7-28

scatter3 is used to plot the data points and this is then overlaid with a plot of the fittedplane using fimplicit3.

请先登录，再进行评论。

Answer 2

Mathieu NOE 2023-7-20

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1998593-how-to-fit-an-equation-catenary-equation-in-2d-to-a-3d-dataset-experimental-data-of-a-catenary-cu#answer_1276408

在 MATLAB Online 中打开

This is a lazy (tired) guy answer

i figure out I could easily transform the 3D array into 2D then do a parabola fit (yes I know it's only an approximation for the catenary cosh model)

  load('export_dataset1.mat')
  x = cat1_2_x';
  y = cat1_2_y';
  z = cat1_2_z';
  m = numel(x);
  figure(1),
  plot3(x,y,z,'*-');
  xlabel('X')
  ylabel('Y')
  zlabel('Z')
  
  P = [x y z]; % create array of x,y,z coordinates
  V = diff(P,1); % array of vectors 
  d = sqrt(sum(V.^2,2)); % array of vectors length
  
  % check distance between first and last point (debug)
  dmax = P(m,:)-P(1,:);
  dmax3D = sqrt(sum(dmax.^2))
  % compute raltive angle between successive vectors
    ang(1) = 0;
    for k = 1:m-2
       ang(k+1) = acos(dot(V(k,:)/norm(V(k,:)),V(k+1,:)/norm(V(k+1,:))));
    end
  
    angc = cumsum(ang);
    angc = angc - angc(end)/2; % so the 2plot shows a balanced (symetrical) start and end angles
    
% create the 2D points by adding vector coordinates
x(1)= 0;
y(1) = 0;
  for k = 1:m-1
    x(k+1) = x(k)+d(k)*cos(angc(k));
    y(k+1) = y(k)+d(k)*sin(angc(k));
  end
  
% check distance between first and last point (debug)
P = [x y];
dmax = P(m,:)-P(1,:);
dmax2D = sqrt(sum(dmax.^2))
% this is the same as the original data
% now you can fit your 2D data
% as I am lazy and it's getting late here I am doing a parabola fit (and
% not a true cosh model fit
%Set up the appropriate matrix A to find the best-fit parabola of the form y=C+Dx+Ex^2.  The 
%first column of A will contain all 1's, using the ones() command.  The second column of A 
%contains x values that are stored in X.  The third column of A contains the squared x values
%that are stored in X.  Elementwise multiplication of X by itself, using .* operator, will 
%produce the desired values for the third column.  
A = [ones(m,1) x x.*x];
A_transposeA = A.' * A;
A_transposeY = A.' * y;
%backslash operation to solve the overdetermined system.
Soln2 = A_transposeA\A_transposeY;
%x values to use for plotting the best-fit parabola.  
xfit = linspace(min(x),max(x),100);
yfit = Soln2(1) + Soln2(2)*xfit + Soln2(3)*xfit.^2; %
figure(2),
plot(xfit, yfit, x, y, 'r*','Markersize',15);
legend('fit','data');
xlabel('X')
ylabel('Y')

3 个评论
显示 1更早的评论隐藏 1更早的评论

Mathieu NOE 2023-7-20

在 MATLAB Online 中打开

this is the same idea for the second data set

I had to do some more work to "isolate" one catenary and apply the same logic

tomorrow I will try to finish my work for the second catenary (shoudn't be too complicated !!)

 
  load('export_dataset2.mat')
  x0 = cat6_1_x';
  y0 = cat6_1_y';
  z0 = cat6_1_z';
  
  figure(1),
  plot3(x0,y0,z0,'*-');
  xlabel('X')
  ylabel('Y')
  zlabel('Z')  
  hold on 
  
  % isolate one catenary data (the one which is parrallel to the Y
  % direction)
   n= find(x0>-0.4 & x0<-0.3);
  x = x0(n);
  y = y0(n);
  z = z0(n);
  
  % sort / unique in y ascending order
  [y,ind,~] = unique(y);
  x = x(ind);
  z = z(ind);
  m = numel(x);
  plot3(x,y,z,'*r','Markersize',15);
  xlabel('X')
  ylabel('Y')
  zlabel('Z')
  hold off
  
  P = [x y z]; % create array of x,y,z coordinates
  V = diff(P,1); % array of vectors 
  d = sqrt(sum(V.^2,2)); % array of vectors length
  
  % check distance between first and last point (debug)
  dmax = P(m,:)-P(1,:);
  dmax3D = sqrt(sum(dmax.^2))
  % compute raltive angle between successive vectors
    ang(1) = 0;
    for k = 1:m-2
       ang(k+1) = acos(dot(V(k,:)/norm(V(k,:)),V(k+1,:)/norm(V(k+1,:))));
    end
  
    angc = cumsum(ang);
    angc = angc - angc(end)/2; % so the 2plot shows a balanced (symetrical) start and end angles
    
% create the 2D points by adding vector coordinates
x(1)= 0;
y(1) = 0;
  for k = 1:m-1
    x(k+1) = x(k)+d(k)*cos(angc(k));
    y(k+1) = y(k)+d(k)*sin(angc(k));
  end
  
% check distance between first and last point (debug)
P = [x y];
dmax = P(m,:)-P(1,:);
dmax2D = sqrt(sum(dmax.^2))
% this is the same as the original data
% now you can fit your 2D data
% as I am lazy and it's getting late here I am doing a parabola fit (and
% not a true cosh model fit
%Set up the appropriate matrix A to find the best-fit parabola of the form y=C+Dx+Ex^2.  The 
%first column of A will contain all 1's, using the ones() command.  The second column of A 
%contains x values that are stored in X.  The third column of A contains the squared x values
%that are stored in X.  Elementwise multiplication of X by itself, using .* operator, will 
%produce the desired values for the third column.  
A = [ones(m,1) x x.*x];
A_transposeA = A.' * A;
A_transposeY = A.' * y;
%backslash operation to solve the overdetermined system.
Soln2 = A_transposeA\A_transposeY;
%x values to use for plotting the best-fit parabola.  
xfit = linspace(min(x),max(x),100);
yfit = Soln2(1) + Soln2(2)*xfit + Soln2(3)*xfit.^2; %
figure(2),
plot(xfit, yfit, x, y, 'r*','Markersize',15);
legend('fit','data');
xlabel('X')
ylabel('Y')
  

Mathieu NOE 2023-7-20

在 MATLAB Online 中打开

well

that was much faster than I thougt

second catenary results

  load('export_dataset2.mat')
  x0 = cat6_1_x';
  y0 = cat6_1_y';
  z0 = cat6_1_z';
  
  figure(1),
  plot3(x0,y0,z0,'*-');
  xlabel('X')
  ylabel('Y')
  zlabel('Z')  
  hold on 
  
  % isolate one catenary data (the one which is parrallel to the Y
  % direction)
   n= find(y0>-0.5 & y0<-0.4);
  x = x0(n);
  y = y0(n);
  z = z0(n);
  
  % sort / unique in x ascending order
  [x,ind,~] = unique(x);
  y = y(ind);
  z = z(ind);
  m = numel(x);
  plot3(x,y,z,'*r','Markersize',15);
  xlabel('X')
  ylabel('Y')
  zlabel('Z')
  hold off
  
  P = [x y z]; % create array of x,y,z coordinates
  V = diff(P,1); % array of vectors 
  d = sqrt(sum(V.^2,2)); % array of vectors length
  
  % check distance between first and last point (debug)
  dmax = P(m,:)-P(1,:);
  dmax3D = sqrt(sum(dmax.^2))
  % compute raltive angle between successive vectors
    ang(1) = 0;
    for k = 1:m-2
       ang(k+1) = acos(dot(V(k,:)/norm(V(k,:)),V(k+1,:)/norm(V(k+1,:))));
    end
  
    angc = cumsum(ang);
    angc = angc - angc(end)/2; % so the 2plot shows a balanced (symetrical) start and end angles
    
% create the 2D points by adding vector coordinates
x(1)= 0;
y(1) = 0;
  for k = 1:m-1
    x(k+1) = x(k)+d(k)*cos(angc(k));
    y(k+1) = y(k)+d(k)*sin(angc(k));
  end
  
% check distance between first and last point (debug)
P = [x y];
dmax = P(m,:)-P(1,:);
dmax2D = sqrt(sum(dmax.^2))
% this is the same as the original data
% now you can fit your 2D data
% as I am lazy and it's getting late here I am doing a parabola fit (and
% not a true cosh model fit
%Set up the appropriate matrix A to find the best-fit parabola of the form y=C+Dx+Ex^2.  The 
%first column of A will contain all 1's, using the ones() command.  The second column of A 
%contains x values that are stored in X.  The third column of A contains the squared x values
%that are stored in X.  Elementwise multiplication of X by itself, using .* operator, will 
%produce the desired values for the third column.  
A = [ones(m,1) x x.*x];
A_transposeA = A.' * A;
A_transposeY = A.' * y;
%backslash operation to solve the overdetermined system.
Soln2 = A_transposeA\A_transposeY;
%x values to use for plotting the best-fit parabola.  
xfit = linspace(min(x),max(x),100);
yfit = Soln2(1) + Soln2(2)*xfit + Soln2(3)*xfit.^2; %
figure(2),
plot(xfit, yfit, x, y, 'r*','Markersize',15);
legend('fit','data');
xlabel('X')
ylabel('Y')
  

VIGNESH BALAJI 2023-7-20

@Mathieu NOE If I am correct you first found a plane where the 2D data lie in the 3D plane and then fitted the 2D data.

If I am correct these lines of code below is same as plane fit -

P = [x y z]; % create array of x,y,z coordinates

V = diff(P,1); % array of vectors

d = sqrt(sum(V.^2,2)); % array of vectors length

% check distance between first and last point (debug)

dmax = P(m,:)-P(1,:);

dmax3D = sqrt(sum(dmax.^2))

% compute raltive angle between successive vectors

ang(1) = 0;

for k = 1:m-2

ang(k+1) = acos(dot(V(k,:)/norm(V(k,:)),V(k+1,:)/norm(V(k+1,:))));

end

angc = cumsum(ang);

angc = angc - angc(end)/2; % so the 2plot shows a balanced (symetrical) start and end angles

% create the 2D points by adding vector coordinates

x(1)= 0;

y(1) = 0;

for k = 1:m-1

x(k+1) = x(k)+d(k)*cos(angc(k));

y(k+1) = y(k)+d(k)*sin(angc(k));

end

请先登录，再进行评论。

Answer 3

John D'Errico 2023-7-20

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1998593-how-to-fit-an-equation-catenary-equation-in-2d-to-a-3d-dataset-experimental-data-of-a-catenary-cu#answer_1276283

编辑：John D'Errico 2023-7-20

How is this a 3-d probem? The curve will lie in a plane. So it is just a 2-d problem. That you have measurements in x,y,z is irrelevant. So reduce it to TWO dimensions first, along the line in the (x,y) plane between the two endpoints.
Then fit the curve, as a function of the new variable along that line between the endpoints. WTP?
Finally, Since you have the equation of the line between the endpoints, convert it back into 3-d, if that is what you really wanted.

3 个评论
显示 1更早的评论隐藏 1更早的评论

Mathieu NOE 2023-7-20

it would be more efficient to share the data

has you tried to find (fit) the plane to which the data belong ?

VIGNESH BALAJI 2023-7-20

Sure, I am sharing the datasets. I am also adding an image of the datasets, you can load the mat files and get their x,y and z co-ordinates.

Dataset 1 conatins one catenary

Dataset 2 contains 2 catenaries

I haven't tried fitting a plane. I was just thinking about it. What could be a good tool to fit a plane, You mean PlanarFit or some other tool to fit a plane ?

请先登录，再进行评论。

Answer 4

VIGNESH BALAJI 2023-7-20

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1998593-how-to-fit-an-equation-catenary-equation-in-2d-to-a-3d-dataset-experimental-data-of-a-catenary-cu#answer_1276428

编辑：VIGNESH BALAJI 2023-7-20

在 MATLAB Online 中打开

@Mathieu NOE thanks a lot for trying this problem even when you are tired.

I also wrote a script in the meantime where I seperated each catenary curve and fitted and fitted a 2d curve,

then i manually rotated and level shiffted the plane to fit the 3d data and to show the results in a 3d fit format.

I have put the code below. This process was very manually intensive as I have to manuallty check the fit.

I think your work showed me an inspiration to find a plane for 3d data and to curve fit in the 2d plane, then I can rotate this plane back to 3d original space to show the fit in 3d space.

I have put down my code for catenary fit (attached images) -

%% curve fitting with catenary equation
% data 
xco = cat1_2_y; % data rearranged and not necessary
yco = cat1_2_z;
zco = cat1_2_x;
%length of catenary used in actual setup
l = 1.76;         % measured from actual cable used in experiment
span = (abs(cat1_2_y(:,1)) + abs(cat1_2_y(:,6)));  % span tells how wide the curve(catenary) is opened up
% finding paramegter "a" to find the sag of the curve 
fun = @(a) (- l/2  + a * sinh(span/(2*a))); x0 = 0.15;
a = abs(fzero(fun,x0));        
xm = linspace(xco(:,1), xco(:,6),100); % divide x data into 100 divisions to get a smooth curve
size_xm = size(xm)
ym = zeros(1,size_xm(2));
for k=1:size_xm(2)
    ym(1,k) = a*cosh((xm(1,k)/a)) + (0.4965-a); % catenary equation in 2d, the plus term is a scaling factor as the curve starts from zero height and actual data is at a particular height
end
zm = linspace(zco(:,1), zco(:,6),100); % divide z data into 100 divisions to get a smooth curve
zm = zm + 0.04;   % 4 cm z error correction, mismatch between experiment and analytical plotting
%% rotate this generated 2d catenary curve in 3d space to fit the 3d experimental data
new = zeros(3,size_xm(2));
for i=1:size_xm(2)
   new(:,i) = rotx(-1.9)*roty(0)*rotz(0)*[xm(:,i); ym(:,i); zm(:,i)]; % x rotation, done manually by seeing plots
end
x_new = new(1,:);
y_new = new(2,:);
z_new = new(3,:);
%% plotting results for comparisoon pf data and fit 
plot3(xco(:), yco(:), zco(:), 'o');    % experimental data
hold on
plot3(xm(:), ym(:), zm(:))              % generated catenary before rotation
xlabel('x axis (m)')
ylabel('y axis (m)')
zlabel('z axis (m)')
plot3(x_new(:), y_new(:), z_new(:), 'g') % generated catenary after rotation
legend('Experimental data', 'Analytical expression before x rotation and static error', 'Analytical expression after x rotation and static error')
length = 2* a*sinh(span/(2*a)); % check length of curve is same as length in experiment
%% fit evaluation, do meansquare error to see goodness of the fit 
% data rearranging
exp_data = [xco;yco;zco];
ana_x = [x_new(:,1), x_new(:,9), x_new(:,23), x_new(:,49), x_new(:,92), x_new(:,100)];
ana_y = [y_new(:,1), y_new(:,9), y_new(:,23), y_new(:,49), y_new(:,92), y_new(:,100)];
ana_z = [z_new(:,1), z_new(:,9), z_new(:,23), z_new(:,49), z_new(:,92), z_new(:,100)];
ana_data = [ana_x; ana_y; ana_z];
error_x = immse(xco, ana_x)
error_y = immse(yco, ana_y)
error_z = immse(zco, ana_z)
% E = rmse( exp_data, new, "all");
 %disp(err)

4 个评论
显示 2更早的评论隐藏 2更早的评论

Mathieu NOE 2023-7-21

well , I don't know howyu could make all this "automatic"

myself I had to do a few manual operations especially for the two catenary data set , to extract each one

then today I have no time for the 3D "backwards" projection topic (sorry)

VIGNESH BALAJI 2023-7-21

@Mathieu NOE I think I didn't explain it clearly. I will explain it now.

I have a script to extract catenaries (that is fine).

I want the fitting part of catenaries in 3D space as automatic instead of manual.

The automatic part of fitting could be a tool or a standard one time written script.

请先登录，再进行评论。

How to fit an equation (catenary equation in 2D) to a 3D dataset (experimental data of a catenary curve in 3D) ?

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

回答（4 个）

6 个评论
显示 4更早的评论隐藏 4更早的评论

3 个评论
显示 1更早的评论隐藏 1更早的评论

3 个评论
显示 1更早的评论隐藏 1更早的评论

4 个评论
显示 2更早的评论隐藏 2更早的评论

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

How to fit an equation (catenary equation in 2D) to a 3D dataset (experimental data of a catenary curve in 3D) ?

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

回答（4 个）

6 个评论 显示 4更早的评论隐藏 4更早的评论

3 个评论 显示 1更早的评论隐藏 1更早的评论

3 个评论 显示 1更早的评论隐藏 1更早的评论

4 个评论 显示 2更早的评论隐藏 2更早的评论

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

WeChat

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

6 个评论
显示 4更早的评论隐藏 4更早的评论

3 个评论
显示 1更早的评论隐藏 1更早的评论

3 个评论
显示 1更早的评论隐藏 1更早的评论

4 个评论
显示 2更早的评论隐藏 2更早的评论