# How can the major axis of a 3d object be determined?

4 views (last 30 days)
Carson Purnell on 7 Jun 2022
Answered: Carson Purnell on 15 Jun 2022
I have some rodlike objects in 3d, below is an example. Is there a way to find its longitudinal axis? I'm trying to generate a bundle of such objects at different orientations and need to know the axis to be able to rotate them about that axis (and properly place them at the correct distance)

Carson Purnell on 15 Jun 2022
Seems like no matlab code has this functionality built in, so I ended up doing it myself. Demo code, assuming you already have a shape within 3d array vol:
%generate points from the volume
[X, Y, Z] = meshgrid(1:size(vol,2), 1:size(vol,1), 1:size(vol,3));
points = [X(:) Y(:) Z(:) vol(:)];
nonzero = points(:,4)>0; %threshold to retain points
filtered = points(nonzero,1:3);
%generate an orientation vector from points
cen = mean(filtered,1);
[~,~,V] = svd((filtered-cen),'econ'); %singular value decomposition, econ for speed
d = V(:,1)'; %the slope vector, ready for input to imrotate3
%demo of line
count =0;
for i=-50:.1:50 %show line
count = count+1;
lin(count,:) = d*i+cen;
end
close all
hold on
plot3(filtered(:,1),filtered(:,2),filtered(:,3),'.')
plot3(lin(:,1),lin(:,2),lin(:,3),'o');
axis equal
function is on the FEX, let me know if it works for anyone else.

### More Answers (1)

Kevin Holly on 9 Jun 2022
Edited: Kevin Holly on 9 Jun 2022
You can figure out the orientation with regionprops3 assuming you have a 3D matrix of the model.
Example:
sample = zeros(100,100,100);
sample(5,:,:)=eye(100);
SE = strel("disk", 3);
new_sample=imdilate(sample,SE);
isosurface(new_sample)
axis equal
xlim([0 100])
ylim([0 100])
zlim([0 100])
rp = regionprops3(new_sample,'Orientation');
rp.Orientation
ans = 1×3
-90.0000 45.0000 0
Let me rotate it with imrotate3
figure
new_sample_rotated = imrotate3(new_sample,90,rp.Orientation,'nearest','loose','FillValues',0);
isosurface(new_sample_rotated)
axis equal
xlim([0 100])
ylim([0 100])
zlim([0 100])
figure
new_sample_rotated = imrotate3(new_sample,35,rp.Orientation,'nearest','loose','FillValues',0);
isosurface(new_sample_rotated)
axis equal
xlim([0 100])
ylim([0 100])
zlim([0 100])
Carson Purnell on 13 Jun 2022
I've found no way to convert this transformation into the real axis of the object, even matlab's method in the aerospace toolbox does not generate anything close to the real axis.