Matlabでクォータニオンのテスト

以下の関数が便利
・quaternion(q0,q1,q2,qe)：コンストラクタ
・conj()：共役クォータニオン。単位クォータニオンの場合は、逆クォータニオンと等しくなる
・parts()：クォータニオンから要素を取り出す。[q0,q1,q2,q3]=parts(quat)のように使う。[,q1,q2,q3]=parts(quat)のようにq0を省略したら動かなかったので注意。

clear; close all;

N=100;
for phi=linspace(0,2*pi,N)
%     phi=pi/3; %60deg;
    u=[1/sqrt(3) 1/sqrt(3) 1/sqrt(3)]; %unit vector(rotational axis)
    q=quaternion(cos(phi/2),u(1)*sin(phi/2),u(2)*sin(phi/2),u(3)*sin(phi/2));
    qinv= q.conj;
    v = [1 0.1 0.1]'; % a vector
    vq =quaternion(0,v(1),v(2),v(3));
    vdashq = q*vq*qinv;
    [v0,v1,v2,v3]=parts(vdashq);
    vdash=[v1 v2 v3]';
    norm(vdash)
    
    quiver3(0,0,0,v(1),v(2),v(3));
    axis equal
    %     view([-10 -10 10]);
    xlim([-1 1]);
    ylim([-1 1]);
    zlim([-1 1]);
    hold on;
    text(0,0,sprintf('phi=%f',phi));
    quiver3(0,0,0,vdash(1),vdash(2),vdash(3));
    quiver3(0,0,0,u(1),u(2),u(3));
    legend('v','vdash','u');
    big;
    drawnow;
    hold off;
end