How to calculate angle difference between 3D vectors and conditionally plot differences
1 view (last 30 days)
Show older comments
I have a 3D image of an object with holes in it. The holes appear more as indents when plotting the image - this is acceptable for me. When looking at the surfnorm and quiver3 plots of the image, I noticed the majority of the surface vectors appear to be at the same angle, while the areas with indents, are at varying angles. I have not been able to verify if this is true or not, however I am working on the assumption that it is true.
I have calculated the normal vectors [Nx, Ny, Nz] using surfnorm, and I have attempted to calculate the angles between each and every vector. The results are not really what I was going for though - the plots look similar to a quiver3 plot of the gradients which makes sense in retrospect.
Below is what I have come up with so far.
% Importing data from .xyz
Data = Drill2Edited;
x = Data(:,1);
y = Data(:,2);
z = Data(:,3);
% Defining step size
xv = linspace(min(x), max(x), 100);
yv = linspace(min(y), max(y), 100);
[X,Y] = meshgrid(xv, yv);
Z = griddata(x,y,z,X,Y);
%% Stem Plots
figure()
subplot(2, 1, 1)
stem3(X, Y, Z)
subplot(2, 1, 2)
stem3(x, y, z)
grid on
xlabel('x')
ylabel('y')
zlabel('z')
%% Surface Plot
figure()
surf(X, Y, Z);
grid on
xlabel('x')
ylabel('y')
zlabel('z')
%% Surface Normal Calculations and Vector Plots
figure()
subplot(2,1,1)
surfnorm(X,Y,Z)
[Nx,Ny,Nz] = surfnorm(X,Y,Z);
subplot(2,1,2)
quiver3(X,Y,Z,Nx,Ny,Nz);
xlabel('x')
ylabel('y')
zlabel('z')
%% Angle Calculations
[theta] = zeros(100);
% All three matrices have the same NaN/non-NaN positions, therefore I am assuming the
% positions are transferrable from matrix to matrix and the positions only
% need to be found in one matrix
% Note: I am only referring to the cell positions, not the values stored
% within those cells/positions
% Creates a for loop iterating through each position in the Nx matrix
for k = 1:numel(Nx)
pos = k;
[r3]=[Nx(pos),Ny(pos),Nz(pos)];
% Uses position loop to find neighbors of every position
for j = pos
szx = size(Nx);
[rx,cx] = ind2sub(szx,pos);
neighx(1:8,1:2) = [rx+[-1;0;1;-1;1;-1;0;1] cx+[-1;-1;-1;0;0;1;1;1] ];
neighx = neighx(all(neighx,2) & neighx(:,1) <= szx(1) & neighx(:,2) <= szx(2),:);
idx = (neighx(:,2)-1)*szx(1) + neighx(:,1); %positions of neighbor
end
% Values of selected positions' neighbors
for i =1:1:length(idx)
[r4]=[Nx(idx(i)),Ny(idx(i)),Nz(idx(i))];
% Calculating Theta
ps = idx(i);
c = cross(r3,r4);
c1 = abs(c);
res = real(acos((dot(r3,r4))/sqrt(c1(1)^2+c1(2)^2+c1(3)^2)));
if (~isnan(res) && res>0)
if (res>theta(ps))
theta(ps) = res;
end
end
end
end
figure()
plot3(theta, Y, Z);
plot3(X, Y, theta); % no indents appear which leads me to believe XY isn't the correct plane to plot on
return
My intended result is a plot with the indented areas and the edges (as I expect the edges to also have changes in the vector angle) highlighted in some way. Eventually, I would like to add a threshold for the vector angle difference but that may be a problem for another time.
I'm using MatLab 2021a
The coordinate points are stored in an xyz file ('Drill2Edited') but it is too large to attach to this question. Instead I have attached 'Drill.mat' with the coordinate data. If another format is preferred let me know.
0 Comments
Answers (0)
See Also
Categories
Find more on Lighting, Transparency, and Shading in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!