MATLAB Answers

How can I retrieve x,y coordinates of a sliced plane through a 3D volume using the obliqueslice function?

조회 수: 22(최근 30일)
Astrid 9 May 2020
댓글: Astrid 21 May 2020
I am using the function obliqueslice to get an aribtrary (not orthogonal) sliced plane of my 3D dataset. I want to get the coordinates of the point that used as input in the resulting image slice. It is necessary to know how the slice is build up. The documentation of the function inlcudes a brief description of the orientation of slice in image plane, but I don't understand exactly how the output image is build up.
My question is: what is meant by "The slice corner that lies close to the origin of the last slice (0,0,P), in the volumetric data constitutes the upper-left pixel in the image plane."?
My interpretation of this explanation is that the slice corner with the shortest Euclidean distance to the point x=0,y=0 and z=size(Volume,3) is always in the left upper corner of the 2D image. However, this is not true in all cases. For example (with example code from documentation):
load mri
V = squeeze(D);
point = [73 5 15.5];
normal = [5 5 -20];
[B,x,y,z] = obliqueslice(V,point,normal);
fig = figure;
axes1 = axes('Parent',fig,'Position',[0.13 0.4 0.37 0.5]);
grid on
view([-62.6 13.8])
title('Position of Slice in 3-D Coordinate Space')
hold on
hold off
legend('Point in the volume','Normal vector','Origin last slice','1','2','3','4','Position',[0.15 0.1 0.3 0.08])
axes2 = axes('Parent',fig,'Position',[0.59 0.39 0.33 0.6]);
title('Output Slice in Image Plane')
corner1 = [x(1,1),y(1,1),z(1,1)];
corner2 = [x(1,end),y(1,end),z(1,end)];
corner3 = [x(end,end),y(end,end),z(end,end)];
corner4 = [x(end,1),y(end,1),z(end,1)];
o_lastslice = [0, 0, size(V,3)];
length_corner1 = sqrt((corner1(1)-o_lastslice(1))^2+(corner1(2)-o_lastslice(2))^2+...
length_corner2 = sqrt((corner2(1)-o_lastslice(1))^2+(corner2(2)-o_lastslice(2))^2+...
length_corner3 = sqrt((corner3(1)-o_lastslice(1))^2+(corner3(2)-o_lastslice(2))^2+...
length_corner4 = sqrt((corner4(1)-o_lastslice(1))^2+(corner4(2)-o_lastslice(2))^2+...
Here, corner3 has the shortest distance to the origin of the last slice, however, corner 1 is taken as upper-left pixel in the image plane.
I am looking forward to get some more information about this function obliqueslice. Other ideas how to get the coordinates of the point in the 2D slice are also welcome. Many thanks!

  댓글 수: 0

Sign in to comment.

채택된 답변

Boubat Matthieu
Boubat Matthieu 11 May 2020
편집: Boubat Matthieu 11 May 2020
I don't know anything about the upper-left pixel of the slice, but this is how I find a point's indices in a slice:
load mri
V = squeeze(D);
pointSubscriptIndicesInVolume = [40 60 13];
normal = [5 5 -20];
[B,x,y,z] = obliqueslice(V,pointSubscriptIndicesInVolume,normal);
pointLinearIndexInSlice = find( round(x) == pointSubscriptIndicesInVolume(1) &...
round(y) == pointSubscriptIndicesInVolume(2) &...
round(z) == pointSubscriptIndicesInVolume(3) );
[pointRow,pointColumn] = ind2sub(size(B),pointLinearIndexInSlice);
This solution requires to use the point's subscript indices instead of its coordinates.
I haven't checked that it will always provide a solution, and only one solution, but it has not failed me so far.

  댓글 수: 1

Astrid 21 May 2020
Hello Boubat Matthieu,
Many thanks for your alternative solution. It works!
I din't always obtain a result in my data, so I'm now using ceil() instead of round() and I'm not only looking at the the pixel equal to the point's subscript indices, but also to one next to it. This sometimes give more than one solution, so I pick the first one (random, but in my case always close enough to the original coordinates):
pointLinearIndexInSlice = find(ceil(x)==ceil(pointSubscriptIndicesInVolume(1)) &...
ceil(y)==ceil(pointSubscriptIndicesInVolume(2)) &...
ceil(z)==ceil(pointSubscriptIndicesInVolume(3)) );
if isempty(pointLinearIndexInSlice)
pointLinearIndexInSlice = find((ceil(x)==ceil(pointSubscriptIndicesInVolume(1))|ceil(x)==(ceil(pointSubscriptIndicesInVolume(1))-1)) &...
(ceil(y)==ceil(pointSubscriptIndicesInVolume(2))|ceil(y)==(ceil(pointSubscriptIndicesInVolume(2))-1)) &...
for j=1:size(pointLinearIndexInSlice,1)
[pointRow,pointColumn] = ind2sub(size(B),pointLinearIndexInSlice);
if size(pointRow,2) > 1
pointRow = pointRow(1);
pointColumn = pointColumn(1);

Sign in to comment.

추가 답변(0개)

Translated by