Auto-completion and smoothing of incomplete 3D ring

hello, here I got an incomplete 3D ring, and I want to interpolate it to a complete and smoothing 3D ring.

Could anyone help me, thanks very much!

Answer 1

[Mathieu NOE]

hello

maybe this ?

load('3D ring.mat')
x = data.Points(:,1);
y = data.Points(:,2);
z = data.Points(:,3);

%% method 1
centroid_x = mean(x);
centroid_y = mean(y);
centroid_z = mean(z);
[theta,r,zz] = cart2pol(x-centroid_x,y-centroid_y,z-centroid_z);


% sort theta in ascending order
[theta,ind] = sort(theta);
r = r(ind);
zz = zz(ind);

% remove duplicates 
[theta,IA,IC] = unique(theta);
r = r(IA);
zz = zz(IA);

theta_new = linspace(min(theta),max(theta),100);
r_new = interp1(theta,r,theta_new);
z_new = interp1(theta,zz,theta_new);

% smoothing
% add some overlap (pre and post data at closure point)
N = 10;
r_tmp = [r_new(end-N:end) r_new r_new(1:1+N)];
rs = smoothdata(r_tmp,'gaussian',N);
zz_tmp = [z_new(end-N:end) z_new z_new(1:1+N)];
zs = smoothdata(zz_tmp,'gaussian',N);
% remove extra pre and post data
rs(1:1+N) = [];
rs(end-N:end) = [];
zs(1:1+N) = [];
zs(end-N:end) = [];

% convert to cartesian
[xn,yn,zn] = pol2cart(theta_new,rs,zs);

% add back centroid info
xn = xn + centroid_x;
yn = yn + centroid_y;
zn = zn + centroid_z;


%% XY plot
figure(1),
plot3(x,y,z,'g.');
hold on
plot3(xn,yn,zn,'r','linewidth',5);
hold off

legend('raw','smoothed');

[Mathieu NOE]

well , Ithink this new code is better to my previous work

In my first submision (above) , if you would remove the smoothing operation, the red curve would show quite a lot of zig zags , which , I admit , is not very good

instead , I prefer this new code - basically the red curve points are generated by another method (extract of the full code) :

% "interpolated/mean" of r,z values : take the average of r and z within an angle domain
Npoints = 100;
theta_new = linspace(min(theta),max(theta),Npoints);
dt = mean(diff(theta_new));

for k = 1:Npoints
    ind = (theta>=theta_new(k)-dt/2) & (theta<theta_new(k)+dt/2);
    if isempty(ind)
        r_new(k) = NaN;
        z_new(k) = NaN;        
    else
        r_new(k) = mean(r(ind));
        z_new(k) = mean(zz(ind));
    end
end

now you get this (a much better result even before we apply any smoothing)

but of course we want to remove nans values and get a closed curve so this is done by adding those 2 lines :

% replace / interpolate the NaN values
r_new = fillmissing(r_new,'linear');
z_new = fillmissing(z_new,'linear');

and we get an already good looking closed curve

now you can add your smoothing on top with smoothn (or any alternative you prefer)

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

% smoothing (if needed)
Y = {xn,yn,zn};
[zz,s,exitflag] = smoothn(Y,1); 
xn = zz{1};
yn = zz{2};
zn = zz{3};

and you get this :

full code :

load('3D ring.mat')
x = data.Points(:,1);
y = data.Points(:,2);
z = data.Points(:,3);
%% method 1
centroid_x = mean(x);
centroid_y = mean(y);
centroid_z = mean(z);
[theta,r,zz] = cart2pol(x-centroid_x,y-centroid_y,z-centroid_z);
% sort theta in ascending order
[theta,ind] = sort(theta);
r = r(ind);
zz = zz(ind);
% remove duplicates 
[theta,IA,IC] = unique(theta);
r = r(IA);
zz = zz(IA);

% "interpolated/mean" of r,z values : take the average of r and z within an angle domain
Npoints = 100;
theta_new = linspace(min(theta),max(theta),Npoints);
dt = mean(diff(theta_new));

for k = 1:Npoints
    ind = (theta>=theta_new(k)-dt/2) & (theta<theta_new(k)+dt/2);
    if isempty(ind)
        r_new(k) = NaN;
        z_new(k) = NaN;        
    else
        r_new(k) = mean(r(ind));
        z_new(k) = mean(zz(ind));
    end
end

% replace / interpolate the NaN values
r_new = fillmissing(r_new,'linear');
z_new = fillmissing(z_new,'linear');

% convert to cartesian
[xn,yn,zn] = pol2cart(theta_new,r_new,z_new);
% add back centroid info
xn = xn + centroid_x;
yn = yn + centroid_y;
zn = zn + centroid_z;

% smoothing (if needed)
Y = {xn,yn,zn};
[zz,s,exitflag] = smoothn(Y,1); 
xn = zz{1};
yn = zz{2};
zn = zz{3};

% closing the curve
xn(end+1) = xn(1);
yn(end+1) = yn(1);
zn(end+1) = zn(1);


%% XY plot
figure(1),
scatter3(x,y,z,'g','MarkerFaceAlpha',.2,'MarkerEdgeAlpha',.2);
hold on
plot3(xn,yn,zn,'r','linewidth',5);
hold off

[Tingting]

Thanks for your patient answer! It really helps me.

[Mathieu NOE]

My pleasure!!