# Need to make some kind of coordinate transformation to 3D cylinder /transformation matrices!

Discussion in 'MATLAB' started by someone, Oct 23, 2011.

1. ### someoneGuest

Hi all

I stole something from this function:

CYLINDER Generate cylinder.
[X,Y,Z] = CYLINDER(R,N) forms the unit cylinder based on the
generator
curve in the vector R.

And now I need help with 3D transformation matrices. Please see the
code below (which will run, if you copy/paste it):
---------------------------------------------------------
function testDELETEme
clc
clear all
close all

figure
axis equal
hold on
grid on
xlabel('x'); ylabel('y'); zlabel('z');
view([-20 70]);

bodies=2; % draw 2 cylinders
hgt=2; % height
bodyPosition=[0 2 -1 ; 1 3 0]';
colors=[rand(3,bodies)];
transformationMatrix{1}=eye(3);
transformationMatrix{2}=rand(3);

for b=1:bodies

% createCylinder(radius, height, xyz_offset, rotation);
try
[X{b}, Y{b}, Z{b}] = createCylinder( bodyRadius(b), ...
hgt, bodyPosition,b), transformationMatrix{b} );
catch
keyboard
end

% e.g. size(sz)=2x21 : 21 x/y-coordinates, for 2 z-values = top
+bottom
sz = size(X{b});

h = surf(X{b}, Y{b}, Z{b});
set(h,'FaceColor',colors,b),'FaceAlpha',0.5);
end
keyboard
end

function [x y z] = createCylinder(radius, height, xyz_offset,
varargin)

if nargin == 4
transformMatrix = varargin{1};
end

n = 20; % The cylinder has "n" points around the circumference.

% Vector "r" contains the radius at equally spaced
% points along the unit height of the cylinder.
r = [radius radius]'; % r = [1 1]';

m = length(r);
if m==1,
r = [r;r];
m = 2;
end
theta = (0:n)/n*2*pi;
sintheta = sin(theta);
sintheta(n+1) = 0;

%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
% NEED TO USE transformMatrix ON THE BELOW x y and z VARIABLES
% IN ORDER TO ROTATE THE CYLINDER!!!
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

% originally, x and y equations didn't have xyz_offset-term in them.
% originally, z = (0:m-1)'/(m-1) * ones(1,n+1);
x = r * cos(theta) + xyz_offset(1);
y = r * sintheta + xyz_offset(2);
z = ((0:m-1)'/(m-1)*height + xyz_offset(3)-0.5*height) * ones(1,n+1);

%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
% NEED TO USE transformMatrix ON THE ABOVE x y and z VARIABLES
% IN ORDER TO ROTATE THE CYLINDER!!!
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

end

---------------------------------------------------------

Problem:

I just cannot figure out what I should do with the 3x3 transformation
matrices, in the end of the code - in order to rotate my 3D
cylinder... Please try the code and you'll see that the cylinders are
not rotated now...

someone, Oct 23, 2011

2. ### someoneGuest

Just to make sure that everyone understand what I'm trying to do: You
know if we have a vector (3x1) we can transform it with a
transformation-matrix (3x3). It goes like this (also see
http://en.wikipedia.org/wiki/Rotation_matrix ):

newVec = transformMat * oldVec

I think you know this... Anyway... I just tried this (see bottom of
the code I copied in here in post #1):

----------------------------------
xvals = r * cos(theta) + xyz_offset(1);
yvals = r * sintheta + xyz_offset(2);
zvals = ((0:m-1)'/(m-1)*height + xyz_offset(3)-0.5*height) * ones(1,n
+1);

%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
% NEED TO USE transformMatrix ON THE ABOVE x y and z VARIABLES
% IN ORDER TO ROTATE THE CYLINDER!!!
%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

for row=1:m
for col=1:n+1
vec = transformMatrix * [xvals(row,col); yvals(row,col);
zvals(row,col)];
x(m,col)=vec(1);
y(m,col)=vec(2);
z(m,col)=vec(3);
end
end
----------------------------------

In my "naive" opinion, that should work. But it doesn't... Can anyone
explain why it doesn't work and maybe tell me what I'm doing wrong?

The solution is probably very simple, the minute somebody posts the
right answer here... I just can't see the solution right now...

someone, Oct 23, 2011

3. ### someoneGuest

Stupid me, 3 lines should be:

x(row,col)=vec(1);
y(row,col)=vec(2);
z(row,col)=vec(3);

And then it "works" for the identity matrix. But everything else gets
screwed up so there's a huge problem with my approach. Something I
cannot see...

someone, Oct 23, 2011
4. ### Matt JGuest

Matt J, Oct 23, 2011
5. ### someoneGuest

Thanks, but I think my problem is really really simple (if you just do
it right and understand it, which I don't but I think others do)... I
prefer not to download a lot of extra code for making a simple problem
work... I think we're talking about that I just need to change a few
lines in my code (the last part) and then I think it'll work - the
code from post number 1 should be able to work with a few changes, I
think... Besides that, it would be interesting for me to understand
why it doesn't work now...

Any clever graphics people in here?

someone, Oct 23, 2011
6. ### Matt JGuest

========

Well, one problem that I think I see is that
transformationMatrix{2}=rand(3);
which is not a valid rotation matrix. A rotation matrix has to be orthogonal with determinant=1.

As a footnote, I would also point out that your method of applying the transformationMatrix, using a for loop is highly inefficient. You can do it much more compactly and efficiently as follows

xyz=[xvals),yvals),zvals)]* transformMatrix.';

xyz=reshape(xyz,[size(xvals),3]);

x=xyz,:,1);
y=xyz,:,2);
z=xyz,:,3);

Matt J, Oct 24, 2011
7. ### Matt JGuest

================

It works fine for me.

Matt J, Oct 24, 2011
8. ### someoneGuest

Aaaah, stupid me
I even did that mistake 2-3 weeks ago, but must've forgotten
everything about that... For now I use:

transformationMatrix{2}=[ ...
0.8599 0.1593 0.4850
-0.5090 0.3406 0.7905
-0.0392 -0.9266 0.3739 ];

det(transformationMatrix{1})
det(transformationMatrix{2})

And first determinant is a pure "1" (for the 3x3 identity), the second
determinant is "1.0000". Is there an easy way for me to check for
orthogonality, so I don't make this mistake again in the future? I use
these 3x3 rotation matrices now and then and it would be nice if I
could make a little check, so the code catches when the rotation
matrix is invalid....?

I guess something with dot product or cross product of the columns
would work, as an orthogonality check - wouldn't it? Or is there a
better way?
I'm always so poor a doing this optimization. Last time I was
interested in this, I head/read that Mathworks were so proud of their
GIT (???? or something ???) compiler, claiming that for loops was no
longer an issue in relation to speed as it had used to be before/
earlier... I'm not calling this routine so much right now (but if I'm
going to make some animation, it's important that I optimize this part
of the code)...

Thanks a lot for your comments... So the only error I did in the end,
was maybe that I didn't came in with a correct rotation matrix... That
was/is ofcourse also a serious thing, but I'm glad that it wasn't any
bigger issue...

someone, Oct 24, 2011
9. ### Matt JGuest

============

You can check if Matrix*Matrix.'=eye(3) to within some tolerance.

================

The JIT does indeed improve performance, but not to the point where it's superior to vectorization. Moreover, the vectorized syntax is shorter and easier to code, so why not prefer it?

Matt J, Oct 24, 2011
10. ### someoneGuest

Great, thanks. I'll try it a bit later!
Ok, thank you. I also prefer vectorization, if I can do it without
thinking too much. Unfortunately, I typically use for-loops because I'm
more used to them.

But thanks a lot for you suggestions/comments. They really helped!

someone, Oct 24, 2011
11. ### someoneGuest

You're right. Works fantastic, thank you

someone, Oct 24, 2011