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

11 ビュー (過去 30 日間)
Carson Purnell
Carson Purnell 2022 年 6 月 7 日
回答済み: Carson Purnell 2022 年 6 月 15 日
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
Carson Purnell 2022 年 6 月 15 日
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.

その他の回答 (1 件)

Kevin Holly
Kevin Holly 2022 年 6 月 9 日
編集済み: Kevin Holly 2022 年 6 月 9 日
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])
  3 件のコメント
Kevin Holly
Kevin Holly 2022 年 6 月 10 日
編集済み: Kevin Holly 2022 年 6 月 10 日
"Euler angles, returned as a 1-by-3 vector. The angles are based on the right-hand rule. regionprops3 interprets the angles by looking at the origin along the x-, y-, and z-axis representing roll, pitch, and yaw respectively. A positive angle represents a rotation in the counterclockwise direction. Rotation operations are not commutative so they must be applied in the correct order to have the intended effect."
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','Centroid');
rp.Orientation
ans = 1×3
-90.0000 45.0000 0
rp.Centroid
ans = 1×3
50.5000 5.0000 50.5000
figure
% new_sample=imtranslate(new_sample,-rp.Centroid);
new_sample_rotated = imrotate3(new_sample,-45,[0 1 0]);
isosurface(new_sample_rotated)
axis equal
xlim([0 100])
ylim([0 100])
zlim([0 100])
rp = regionprops3(new_sample_rotated,'Orientation','Centroid');
rp.Orientation
ans = 1×3
90 0 0
rp.Centroid
ans = 1×3
71.5000 5.0000 71.5000
figure
% new_sample_rotated=imtranslate(new_sample,-rp.Centroid);
new_sample_rotated = imrotate3(new_sample_rotated,90,[0 0 1]);
isosurface(new_sample_rotated)
axis equal
xlim([0 100])
ylim([0 100])
zlim([0 100])
rp = regionprops3(new_sample_rotated,'Orientation');
rp.Orientation
ans = 1×3
0 0 0
Using rotate:
figure
isosurface(new_sample_rotated)
axis equal
xlim([0 100])
ylim([0 100])
zlim([0 100])
h=gca;
h.Children(2)
ans =
Patch with properties: FaceColor: 'flat' FaceAlpha: 1 EdgeColor: 'none' LineStyle: '-' Faces: [5148×3 double] Vertices: [2576×3 double] Show all properties
rotate(h.Children(2),[1 0 0],45)
figure
isosurface(new_sample_rotated)
hold on
isosurface(new_sample_rotated)
axis equal
xlim([0 100])
ylim([0 100])
zlim([0 100])
h=gca;
rotate(h.Children(1),[1 0 0],35)
Carson Purnell
Carson Purnell 2022 年 6 月 13 日
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.

サインインしてコメントする。

カテゴリ

Help Center および File ExchangeImage Segmentation and Analysis についてさらに検索

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by