Problem getting slices of a volume data set along an arbitrary orientation
8 views (last 30 days)
Show older comments
I have a 3D array of uint8 values, and a direction vector. I want to be able to pick out X slices along this vector, with the orientation being the direction the slicer is looking (that is to say, the X and Y axes of the slice are two normals to the direction vector that are perpendicular to each other).
I have the following solution, but the results show at a 90-degree offset to what is expected along the Y axis local to the direction vector, and I cannot work out why.
%Get "up" and "right" direction vectors
up = get_rotated_direction(direction, 90,0,0);
right = get_rotated_direction(direction, 0,90,0);
z = 0;
for i=1:zStep:zCount
z = z +1;
%Slicer follows a path denoted by "points"
point = points(i,:);
%topLeft holds the coordinates in the data volume of the point at [X=0,Y=0] in the slice.
topLeft = point - right * width / 2 - up * height / 2;
for x=1:width
for y=1:height
offsetPoint = round(topLeft + up*y+ right*x);
result(x,y,z) = data(offsetPoint(1),offsetPoint(2),offsetPoint(3));
end
end
end
The result is a valid image, but at the wrong orientation, so I believe the issue lies in get_rotated_direction or my calls to it, however in testing with the vectors [0,0,1], [0,1,0] and [1,0,0] the results appear to be correct.
get_rotated_direction looks like this:
function result = get_rotated_direction (direction, x_angle,y_angle,z_angle)
%get the rotation from origin to direction
rot = vrrotvec([0,0,1],direction);
%Turn this into a quaternion
rotquat = angle2quat(rot(1)*rot(4),rot(2)*rot(4),rot(3)*rot(4),'XYZ');
%get the angles in radians
x_angle = deg2rad(x_angle);
y_angle = deg2rad(y_angle);
z_angle = deg2rad(z_angle);
%Create the desired rotation from the origin
desiredRot = angle2quat(x_angle,y_angle,z_angle,'XYZ');
%Rotate the origin to our desired rotation
%Then rotate that by the rotation from origin to our direction
%This should produce a direction that has been locally rotated.
direction= quatrotate(rotquat,quatrotate(desiredRot,[0,0,1]));
%Normalize the direction and return
result = direction/norm(direction);
end
0 Comments
Answers (1)
Sean de Wolski
on 29 Dec 2014
You should be able to do this by applying the operations to a meshgrid of values. Generate a grid for your coordinates and then generate the z grid as a function of the x/y grid. How you get z from x/y will be based on your orientation+offset.
For example:
% Data
[x, y, z] = meshgrid(-3:0.25:3);
v = x.*exp(-x.^2 - y.^2 - z.^2);
Slice along Z = -X + Y
% Data defining a surface
[xs, ys] = meshgrid(-3:0.25:3);
zs = -xs + ys;
% Slice along it
slice(x, y, z, v, xs, ys, zs)
colorbar
See Also
Categories
Find more on Red 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!