Finding angle of rotations from a given unit vector to rotate a given vector using those angle to align with the previous unit vector

I have a unit vector , let's say, b=[0.1844 -0.7417 0.6449] and I want to find the angle of rotations so that I can align a given vector a=[0 0 1] onto that unit vector,b. The way I started this problem is by finding the angle of rotations from the vector and then plugging them back inti the rotation matrix to allign a onto b in the same direction. But I'm stuck at this point. It's not giving me the angles correctly. I've attached my code below.

a=[0,0,1];

[R,theta_x,theta_y,Ry,Rx]=rot_xy(a(1,1),a(1,2),a(1,3));

test=R*a';
function [w,w_k,w_l,w_jk,w_kl ] = rot_xy( x,y,z)
% theta_x=-atand(y/sqrt(x^2+z^2));
% theta_y=pi-atand(x/z);

w_k=-atan(y/sqrt(x^2+z^2));
w_l=pi-atan(x/z);

%rot about y: tilt rotation
w_jk = [  cos(w_k)     0      sin(w_k) ;
             0         1         0     ;
         -sin(w_k)     0      cos(w_k) ];

%rot about x: twist rotation 
w_kl =  [     1         0         0      ;
              0      cos(w_l) -sin(w_l)  ;
              0      sin(w_l)  cos(w_l)  ]; 
        
         
   w = w_jk * w_kl;
end
figure
 plot3([0,a(1,1)],[0,a(1,2)],[0,a(1,3)],'r')
 grid on
 ax=gca;
 ax.XColor='green';
 ax.YColor='magenta';
 grid on
 axis([-.5 .5 -2.5 1 -1.5 1.5])
 hold on
 plot3([0,1.5*(b(1,1))],[0,1.5*(b(1,2))],[0,1.5*(b(1,3))],'k')
 hold on
 c=test';
hold on
 plot3([0,c(1,1)],[0,c(1,2)],[0,c(1,3)],':r')

[Matt J]

Generally speaking, representing the alignment as a rotation about x followed by y will either be impossible, or have two solutions. For example, it is easy to see that there is no solution for

a=[1 0 0 ].';
b=[0 1 0].';

Direct substitution into the equation b = w_jk * w_kl *a leads to

 [0 1 0] = [cos(w_k)     0      -sin(w_k)]

which clearly has no solution.

Conversely, the following data has two solutions and two correspondingly different rotation matrices. You can reach b from a with either a rotation about x of t degrees followed by a 78 degree rotation about y. Or, you can do a rotation about x of -t degrees followed by a rotation about y of 102 degrees.

 a =[

  -0.874881810907645
   0.484336470795830
   0.000000000000000]


b =[

                   0
   0.447213595499958
   0.894427190999916]

t =  22.578683389390513;

Answer 1

[Matt J]

The following might be what you're looking for. It uses my AxelRot utility from the File Exchange

https://www.mathworks.com/matlabcentral/fileexchange/30864-3d-rotation-about-shifted-axis?focused=5191309&tab=function

Or, if you only want theta, you could modify the code to return only that, which then wouldn't require AxelRot.

function [M,theta, Nrm]=vecrot(vstart,vend)
%Find rotation carrying one vector toward another about their common perpendicular
%axis.
%           
%IN:         
%            
%    vstart: Initial vector
%      vend: Final vector
%            
%OUT:        
%            
%    M:  homogeneous 4x4 rotation matrix carrying vstart to vend in a
%        rotation about the axis Nrm=cross(vstart,vend)
%    theta:  the rotation angle in degrees
%    Nrm: the rotation axis



  vstart=vstart(:)/norm(vstart);
  vend=vend(:)/norm(vend);
 
  Nrm=cross(vstart,vend);
   
  b=vend.'*vstart;
  
  theta = atan2d(sqrt(1-b^2),b);

   M=AxelRot(theta,Nrm,[]);

[Nafisa Amin]

But it's giving me only one angle of rotation. But I need two rotational angle, one for y-axis and one for x-axis. How to find both the angles?

[Matt J]

Why do you need 2 angles?

[Nafisa Amin]

One for rotating in y axis and another for x axis so that I can use those in my rotation matrix. I need the tilt and twist angle from my uploaded code

[Matt J]

But the code I gave you gives you the rotation matrix already

a=[0 0 1]
b=[0.1844 -0.7417 0.6449];


[M,theta, Nrm]=vecrot(a,b);
R=M(1:3,1:3);

%%test

>> R*a(:),b(:)

ans =

    0.1844
   -0.7417
    0.6449


ans =

    0.1844
   -0.7417
    0.6449

And regardless, how do you know that it is possible to obtain the alignment in the two-angle form you have posted? I gave an example in my comment which shows that it may be impossible for certain a and b.

[Nafisa Amin]

Thanks a lot for your help. I was wondering if this Axelrot utility can be used for rotating a matrix. For my research project, I need to rotate an initial matrix to another matrix and I was trying to use this utility but it's giving me error. suppose, my initial matrix is, E=[1,0,0;0,2,0;0,0,4] and I want this to end in E_1 = [0.25, 0, 0; 0, 0.25, 0; 0, 0, .497];

This is what I've done so far,

e=E(:)/norm(E);
e1=E_1/norm(E_1);
nrm=cross(e,e1);
theta = atan2d(sqrt(1-bb^2),bb);
M=AxelRot(theta,nrm,[])

%%%%%%%%%%%%%%%%%%%%%%%%

Error using  * 
Incorrect dimensions for matrix
multiplication. Check that the number
of columns in the first matrix matches
the number of rows in the second
matrix. To perform elementwise
multiplication, use '.*'.

Error in AxelRot (line 85)
    AxisShift=x0-(x0.'*u).*u;

%%%%%%%%%%%%%%%%%%%%%%

I tried to fix it but still giving me error,

Matrix dimensions must agree.

Error in AxelRot (line 85)
    AxisShift=x0-(x0.'.*u).*u;

%%%%%%%%%%%%%%%%%%%%%%

I would really appreciate your help on this regard.

[Matt J]

I do not see in what way E_1 is a rotation of E.

[Nafisa Amin]

I was thinking of an arbitrary matrix for the initial one. But my final rotated vector should be E_1. Is there any way to do that?

[Matt J]

I do not understand what you mean by "rotating" a matrix unless you mean what is done by imrotate

https://www.mathworks.com/help/images/ref/imrotate.html

[Nafisa Amin]

Thanks for all of your help. I really appreciate.

[Nafisa Amin]

Hi, I am facing one problem for which I need help. If I have two unit bectors a and b and if a=b or a=-b, the file you have given is showing NaN.

a=[0 0 1];
b=[0 0 -1];
[M,theta, Nrm]=vec_rot(a,b);
R=M(1:3,1:3);

Warning: Matrix is singular, close to singular or badly scaled.
Results may be inaccurate. RCOND = NaN. 

> In AxelRot (line 94)
  In vec_rot (line 24) 
  
>> rot=R*a(:);
>> rot
rot =

   NaN
   NaN
   NaN

[Matt J]

How would you define the result in those cases? There are infinite choices, but this one may serve:

if all(a==b)
    
    R=eye(3);
    
elseif all(a==-b)
    
   N=null(a(:).');
   c=N(:,1);
   M=vec_rot(a,c);
   R=M(1:3,1:3)^2; 
   
else
    
   [M,theta, Nrm]=vec_rot(a,b);
   R=M(1:3,1:3);
    
end

[Nafisa Amin]

Thanks. It works.