3D curve fitting

I am a beginner in MATLAB, and now I have obtained a point cloud data for 3D curve fitting relative to these points, not surface fitting. Is there any method that can achieve good 3D curve fitting? thanks

[Mathieu NOE]

you could look at this (File exchange)

polyfitn - File Exchange - MATLAB Central (mathworks.com)

[tabf]

Thank you very much for your reply, but the main function I want to achieve now is as follows: curve fitting on point cloud data, and then interpolating the obtained curve. For example, the number of interpolation points is t, and finally using t to represent the expressions of x, y, and z

[Mathieu NOE]

do you have started a code ? do you have some data ?

[tabf]

I have data, which is a set of point clouds with approximately 800 points, and the curve type is S-shaped. I also attempted to write a code that utilized the spline function, but the results were not very good. Moreover, the polynomial coefficients of each data point were returned instead of the equation coefficients of a curve. Do you have any good suggestions. Sincerely thank you for your reply.

[Mathieu NOE]

do you mind sharing your code and data ?

[tabf]

ptCloud = pcread('quxiandian.pcd');

N = ptCloud.Location;

x = N(:, 1);

y = N(:, 2);

z = N(:, 3);

p_xy = polyfit(y, x, 4);

p_xz = polyfit(x, z, 4);

y_new = linspace(min(y), max(y), 500);

x_new = polyval(p_xy, y_new);

z_new = polyval(p_yz, x_new);

xyz_new = [x_new', y_new', z_new'];

figure;

view(3);

line(xyz_new(:, 1), xyz_new(:, 2), xyz_new(:, 3), 'Color','b', 'LineWidth', 2);

grid on;

[Mathieu NOE]

I wonder if you want to fit a model or simply smooth the data and get something like this :

if this is what you want , simply download this FEX submission :

smoothn - File Exchange - MATLAB Central (mathworks.com)

and use this code

x = N(:, 1);
y = N(:, 2);
z = N(:, 3);

  u = smoothn({x,y,z},1e4);
  plot3(x,y,z,'r.',u{1},u{2},u{3},'k','linewidth',2)
  axis tight square

[tabf]

Thank you very much for your reply. I would like to obtain a model to represent these points, which can be segmented and ultimately output the parameters of this model

[tabf]

Thank you for your help, but I found that the existing methods have great differences in the effect of fitting different data. Here are my original data and the functions I want to achieve: we can clearly see that there is an S-shaped groove in the middle of this device. Now the purpose is to propose the bottom characteristic line of the S-shaped groove and perform curve fitting. After obtaining the fitting curve, I think that interpolation is still needed, Because in the end, I want to send this curve to the robot. I hope those who have the ability to understand this knowledge can provide me with some help. Finally, I sincerely thank everyone for their replies once again.

[Mathieu NOE]

this is a code to find a polynomial fit for the S shaped groove (trajectory)

N = readmatrix('S.txt');
x = N(:, 1);
y = N(:, 2);
z = N(:, 3);

% detrend the Z data
order = 1;
p = polyfitn([x,y],z,order);
pC = p.Coefficients; % get the polynomial coefficients
pTerms = p.ModelTerms;

% create the polynomial model (z = f(x,y))
zt = 0;
for k = 1:numel(pC)
    zt = zt + pC(k)*(x.^pTerms(k,1)).*(y.^pTerms(k,2)); %
end

figure(1),
  plot3(x,y,z,'r.',x,y,zt,'.k','linewidth',2);  % 
  xlabel('X');
  ylabel('Y');
  zlabel('Z');
  legend('raw data','fitted plane');
  axis tight square

% apply detrend to the Z data
zd = z - zt;
figure(2),
  plot3(x,y,zd,'.','linewidth',2);  % 
  xlabel('X');
  ylabel('Y');
  zlabel('Z');
  axis tight square


% keep the highets z points to get the S shape of the groove
id = (zd>0.85*max(zd));
xx = x(id);
yy = y(id);
% make sure x data is unique and sorted
[xx,ia,ic] = unique(xx);
yy = yy(ia);

% Fit a polynomial p of degree "degree" to the (x,y) data:
degree = 5;
p = polyfit(xx,yy,degree);

% Evaluate the fitted polynomial p and plot:
yyf = polyval(p,xx);
eqn = poly_equation(flip(p)); % polynomial equation (string)
Rsquared = my_Rsquared_coeff(yy,yyf); % correlation coefficient

figure(3);plot(xx,yy,'*',xx,yyf,'-')
  xlabel('X');
  ylabel('Y');
legend('data',eqn)
title(['Data fit , R² = ' num2str(Rsquared)]);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function Rsquared = my_Rsquared_coeff(data,data_fit)
    % R² correlation coefficient computation
    
    % The total sum of squares
    sum_of_squares = sum((data-mean(data)).^2);
    
    % The sum of squares of residuals, also called the residual sum of squares:
    sum_of_squares_of_residuals = sum((data-data_fit).^2);
    
    % definition of the coefficient of correlation is
    Rsquared = 1 - sum_of_squares_of_residuals/sum_of_squares;
    
    
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function eqn = poly_equation(a_hat)
eqn = " y = "+a_hat(1);
for i = 2:(length(a_hat))
    if sign(a_hat(i))>0
        str = " + ";
    else
        str = " ";
    end
    if i == 2
        eqn = eqn+str+a_hat(i)+" * x";
    else
        eqn = eqn+str+a_hat(i)+" * x^"+(i-1)+" ";
    end
end
eqn = eqn+" ";
end
        

[tabf]

I tried to fit the 2D curves of XY and yz first, and then merge the two curves into one 3D curve, but the results I have obtained are not very good now

Answer 1

[Mathieu NOE]

hello again

I can make you this suggestion

I devised that your curve could be parametrized by these 2 equations :

z = a + b*y;
x = c + d*sin(w*y+e);

hope it helps

% ptCloud = pcread('quxiandian.pcd'); 
% N = ptCloud.Location;

N = readmatrix('QUXIANDIAN.txt');

x = N(:, 1);
y = N(:, 2);
z = N(:, 3);

%% model fit 
% z = a + b*y;
% x = c + d*sin(w*y+e);

% initial values (for further optimisation - see below)
w = (2*pi)/(max(y)-min(y));
b = (max(z)-min(z))/(max(y)-min(y));
a = z(1) - b*y(1);
c = mean(x);
d = (max(x)-min(x))/2;
e = 4;

% create data for the fit model
yf = linspace(min(y),max(y),100);
zf = a + b*yf;
xf = c + d*sin(w*yf+e);


% Fit a polynomial p of degree "degree" to the (y,z) data:
    degree = 1;
    p = polyfit(y,z,degree);
    a = p(2);
    b = p(1);

% Fit custom equation to the (x,y) data:
    % option 1 : with fminsearch
    f = @(c,d,e,w,y) c + d*sin(w*y+e);
    obj_fun = @(params) norm(f(params(1), params(2), params(3), params(4), y)-x);
    C1_guess = [c d e w];
    sol = fminsearch(obj_fun, C1_guess); % 

    % update c,d,e,w
    c = sol(1);
    d = sol(2);
    e = sol(3);
    w = sol(4);
    zf = a + b*yf;
    xf = c + d*sin(w*yf+e);
    figure(1),plot3(x,y,z,'r.',xf,yf,zf,'k','linewidth',2)
    axis tight square

[Mathieu NOE]

hello again

FYI, you can also do a polynomial fit using this excellent FEX submission :

polyfitn - File Exchange - MATLAB Central (mathworks.com)

code :

N = readmatrix('QUXIANDIAN.txt');

x = N(:, 1);
y = N(:, 2);
z = N(:, 3);00;

p = polyfitn([x,y],z,3);
% % FEX : https://fr.mathworks.com/matlabcentral/fileexchange/34765-polyfitn?s_tid=ta_fx_results

% % The result can be converted into a symbolic form to view the model more simply. 
% % Here I'll use the sympoly toolbox, but there is also a polyn2sym function provided.
% % FEX : https://fr.mathworks.com/matlabcentral/fileexchange/9577-symbolic-polynomial-manipulation?s_tid=srchtitle
% polyn2sympoly(p)
% % ans =
% %     0.0011322*X1^3 + 0.0010727*X1^2*X2 - 0.28262*X1^2 - 0.00058434*X1*X2^2 - 0.10892*X1*X2 + 20.7666*X1 - 0.00022656*X2^3 + 0.062697*X2^2 + 2.5926*X2 - 121.0331

p = p.Coefficients; % get the polynomial coefficients

% create clean smooth x,y data 
yf = linspace(min(y),max(y),200);
xf = interp1(y,x,yf);
% smooth a bit xf
xf = smoothdata(xf,'gaussian',10);

% create the polynomial model (z = f(x,y))
zf = p(1)*xf.^3 + p(2)*xf.^2.*yf + p(3)*xf.^2 + p(4)*xf.*yf.^2 + p(5)*xf.*yf + p(6)*xf + p(7)*yf.^3 + p(8)*yf.^2 + p(9)*yf + p(10);

  plot3(x,y,z,'r.',xf,yf,zf,'k','linewidth',2)
  axis tight square