dimanche 10 février 2013

Rotation 3D autour d'un axe


Les fichiers "rotation3d_cube.pdf" et "rotation3d_cube.tex" sont dans le dossier :


En utilisant les quaternions :



Les fichiers "rotation3d-q_2.pdf" et "rotation3d-q_2.pdf" sont dans le dossier indiqué au début.

La macro  \psRotationIIID permet de fondre les 2 précédentes en une par l'introduction d'un booléen [quaternion] dont le positionnement permet de choisir la méthode.

Cette macro \psRotationIIID est un complément de pst-solides3d. Elle utilise comme paramètres :
  • [angle=valeur], qui est l’angle de la rotation en degrés ;
  • [axe=x y z], qui doit contenir les 3 composantes du vecteur directeur de l’axe en coordonnées
    cartésiennes. Si on préfère les coordonnées sphériques, on fera suivre celles-ci de
    rtp2xyz : [r rtp2xyz] ;
  • et un booléen [quaternion] positionné par défaut à false. Par défaut le calcul s’effectue
    avec la matrice de rotation 3D classique. Si le booléen est activé par la simple écriture de
    [quaternion] dans les options, le calcul s’effectue avec la méthode des quaternions. 
Cette courte documentation contient aussi un bref aperçu de la théorie des quaternions et de leur lien avec les rotations en 3D écrit par Michel Petitjean.
Les fichiers "rotation3d/rotationIIID.pdf" et "rotation3d/rotationIIID.pdf" sont dans le dossier :


\documentclass[12pt]{article}
\usepackage[a4paper,margin=2cm]{geometry}
\usepackage[latin1]{inputenc}
\usepackage[T1]{fontenc}
\usepackage[garamond]{mathdesign}
\usepackage{pst-solides3d}
\usepackage[colorlinks=true]{hyperref}
\usepackage{amsmath}
\makeatletter
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% la macro \psRotationIIID
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Rotation autour d'un axe passant par O
\define@key[psset]{pst-solides3d}{angle}{\def\psk@solides@angle{#1 }}% rotation autour de l'axe en degres
\psset[pst-solides3d]{angle=0}
%
\define@boolkey[psset]{pst-solides3d}[Pst@]{quaternion}[true]{}
\psset{quaternion=false}%
%\psRotationIIID(x,y,z)
% les coordonnées du vecteur axe sont données par axe=ux uy uz
\def\psRotationIIID{\pst@object{psRotationIIID}}
\def\psRotationIIID@i(#1,#2,#3){{%
  \pst@killglue%
  \use@par%
\addto@pscode{
    \tx@optionssolides
    SolidesDict begin
    /Angle \psk@solides@angle def
    /cosT {Angle cos} bind def
    /sinT {Angle sin} bind def
\ifPst@quaternion
    /Angle_2 {Angle 2 div} bind def
    /p_ {Angle_2 cos} bind def
    /sinA_2 {Angle_2 sin} bind def
    axe unitaire3d
    /u3 exch sinA_2 mul def /u2 exch sinA_2 mul def /u1 exch sinA_2 mul def
    /t2  {p_ u1 mul} bind def
    /t3  {p_ u2 mul} bind def
    /t4  {p_ u3 mul} bind def
    /t5  {u1 dup mul neg} bind def
    /t6  {u1 u2 mul} bind def
    /t7  {u1 u3 mul} bind def
    /t8  {u2 dup mul neg} bind def
    /t9  {u2 u3 mul} bind def
    /t10 {u3 dup mul neg} bind def
/Rotation3D {
4 dict begin
   /M defpoint3d % on récupère les coordonnées
   M /z exch def
     /y exch def
     /x exch def
 t8 t10 add x mul t6 t4 sub y mul add t3 t7 add z mul add 2 mul x add % x'
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 t4 t6 add x mul t5 t10 add y mul add t9 t2 sub z mul add 2 mul y add % y'
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 t7 t3 sub x mul t2 t9 add y mul add t5 t8 add z mul add 2 mul z add  % z'
end} def
\else
   axe unitaire3d
    /u3 exch def /u2 exch def /u1 exch def
/Rotation3D {
5 dict begin
   /M defpoint3d % on récupère les coordonnées
   M /z exch def
     /y exch def
     /x exch def
     /S u1 x mul u2 y mul add u3 z mul add def
   u1 S mul 1 cosT sub mul x cosT mul add u2 z mul u3 y mul sub sinT mul add % x'
   u2 S mul 1 cosT sub mul y cosT mul add u3 x mul u1 z mul sub sinT mul add % y'
   u3 S mul 1 cosT sub mul z cosT mul add u1 y mul u2 x mul sub sinT mul add % z'
end} def
\fi
    end}%
\psSolid[transform=Rotation3D](#1,#2,#3)
}\ignorespaces%
}
\makeatother
\title{Rotation 3D autour d'un axe passant par l'origine}
\date{Février 2013}
\begin{document}

\maketitle
La macro \verb+\psRotationIIID+ est un complément de \textsf{pst-solides3d}. Elle utilise comme paramètres~:
\begin{itemize}
  \item \texttt{[angle=valeur]}, qui est l'angle de la rotation en degrés ;
  \item \texttt{[axe=x y z]}, qui doit contenir les 3 composantes du vecteur directeur de l'axe en coordonnées cartésiennes. Si on préfère les coordonnées sphériques, on fera suivre celles-ci de \texttt{rtp2xyz}~: [r $\theta$ $\phi$ rtp2xyz] ;
  \item et un booléen \texttt{[quaternion]} positionné par défaut à \texttt{false}. Par défaut le calcul s'effectue avec la matrice de rotation 3D classique. Si le booléen est activé par la simple écriture de \textsf{[quaternion]} dans les options, le calcul s'effectue alors avec la méthode des quaternions.
\end{itemize}
\section{Calculs avec la matrice}
\begin{verbatim}
\psRotationIIID[object=cube,
            name=A,
            a=2,
            axe=0 1 1,
            angle=-60,
            RotX=45,
            fillcolor=white](0,3 2 sqrt mul,3 2 sqrt mul)
\end{verbatim}
\begin{center}
\begin{pspicture}(-7,-1)(7,6)
\psset{viewpoint=50 20 20 rtp2xyz,Decran=50,lightsrc=viewpoint}
\psset{solidmemory}
\psLineIIID[linestyle=dashed](0,0,0)(0,5,5)
\psRotationIIID[object=cube,
            name=A,
            a=2,
            axe=0 1 1,
            angle=-60,
            RotX=45,
            fillcolor=white](0,3 2 sqrt mul,3 2 sqrt mul)
\multido{\i=0+1}{6}{
\psSolid[object=plan,action=none,
         definition=solidface,args=A \i,name=P\i]}
\psProjection[object=texte,linecolor=red,text=D,plan=P0,phi=-90,fontsize=50]%
\psProjection[object=texte,linecolor=red,text=A,plan=P1,phi=90,fontsize=50]%
\psProjection[object=texte,linecolor=red,text=C,plan=P2,phi=-90,fontsize=50]%
\psProjection[object=texte,linecolor=red,text=F,plan=P3,phi=90,fontsize=50]%
\psProjection[object=texte,linecolor=red,text=E,plan=P4,phi=180,fontsize=50]%
\psProjection[object=texte,linecolor=red,text=B,plan=P5,phi=180,fontsize=50]%
\composeSolid
\psSolid[object=vecteur,
        definition={[.1 .5]},
        linecolor={[cmyk]{1,0,1,0.5}},
        args=0 1 1](0,5,5)
\psSolid[object=vecteur,
        definition={[.1 .5]},
        linecolor=red,
        args=0 0 2](0,0,0)
\psSolid[object=vecteur,
        definition={[.1 .5]},
        linecolor=blue,
        args=0 2 0](0,0,0)
\psSolid[object=vecteur,
        definition={[.1 .5]},
        linecolor=green,
        args=2 0 0](0,0,0)
\end{pspicture}
\end{center}
\section{Calculs avec le quaternion}
\begin{verbatim}
\def\axe{5 5 0 }
\psRotationIIID[object=cylindre,h=6,r=0.4,
            name=A,
            fcol=2 1 33 { (rouge) } for,
            axe=\axe,
            angle=-40,
            quaternion,
            fillcolor=white,ngrid=4 32](3,3,-3)
\end{verbatim}

\def\axe{5 5 0 }
\begin{center}
\begin{pspicture}(-8,-4)(8,4)
\psset{viewpoint=50 30 20 rtp2xyz,Decran=50,lightsrc=viewpoint}
\psLineIIID[linestyle=dashed](0,0,0)(3.28,3.28,0)
\psset{solidmemory}
\psSolid[object=plan,
         definition=normalpoint,
         args={0 0 0 [0 0 1 180]},
         action=draw,linecolor=red,
         planmarks,
         base=-7 7 -7 7]
\psSolid[object=vecteur,
        definition={[.1 .5]},
        linecolor=red,
        args=0 0 1](0,0,0)
\psSolid[object=vecteur,
        definition={[.1 .5]},
        linecolor=blue,
        args=0 1 0](0,0,0)
\psSolid[object=vecteur,
        definition={[.1 .5]},
        linecolor=green,
        args=1 0 0](0,0,0)
\psRotationIIID[object=cylindre,h=6,r=0.4,name=A,
            fcol=2 1 33 { (rouge) } for,
            axe=\axe,
            angle=-40,
            quaternion,
            fillcolor=white,ngrid=4 32](3,3,-3)
\psSolid[object=plan,action=none,
   definition=solidface,args=A 0,name=P0]
\psSolid[object=plan,action=none,
   definition=solidface,args=A 1,name=P1]
\psProjection[object=texte,linecolor=red,text=A,plan=P0,phi=90,fontsize=20]%
\psProjection[object=texte,linecolor=red,text=B,plan=P1,phi=90,fontsize=20]%
\composeSolid
\psSolid[object=vecteur,
        definition={[.1 .5]},  %% radius height
        linecolor={[cmyk]{1,0,1,0.5}},
        args=2 2 0](3.28, 3.28, 0)
\psPolygonIIID[linecolor=red](7,7,0)(7,-7,0)(-7,-7,0)(-7,7,0)
\end{pspicture}
\end{center}
\section{Quelques éléments de théorie}
Cet aperçu de la théorie des quaternions et de leur lien avec les rotations en 3D a été écrit par Michel Petitjean\footnote{\url{http://petitjeanmichel.free.fr/itoweb.petitjean.html}}, que je remercie sincèrement pour sa contribution.

Le quaternion unitaire est à la rotation 3D ce qu'un nombre
complexe de module 1 est à la rotation 2D.

Dans le cas général (non unitaire), le quaternion peut s'écrire
$q = p + a\mathbf{i} + b\mathbf{j} + c\mathbf{k}$. $p$, $a$, $b$ et $c$ étant des réels, et $\mathbf{i}$, $\mathbf{j}$,
et $\mathbf{k}$ étant des imaginaires vérifiant : $\mathbf{i}^2 = -1$, $\mathbf{j}^2 = -1$, $\mathbf{k}^2 =
-1$ (pour un complexe on aurait seulement $p + a\mathbf{i}$). Il faut en
plus définir les règles des produits entre ces 3 types
d'imaginaires. Entre un réel et un des imaginaires, c'est comme
d'habitude. Entre les complexes, on a : $\mathbf{ij} = - \mathbf{ji} = \mathbf{k}\quad \mathbf{jk} = -
\mathbf{kj} = \mathbf{i}\quad \mathbf{ki} = - \mathbf{ik} = \mathbf{j}$. Maintenant avec ces règles on
sait définir le produit de 2 quaternions.

%Notation: l'apostrophe désignera une transposition (vecteur ou matrice).

Tout comme un complexe de module 1 représente une rotation 2D,
le quaternion unitaire $(p^2+a^2+b^2+c^2=1)$ représente une
rotation 3D, d'axe $\overrightarrow{u}$, avec $\overrightarrow{u}=\left[\begin{array}[m]{l}
a\\
b\\
c
\end{array}\right]$ et d'angle $\theta$ (défini
dans le plan orthogonal à l'axe et orienté par l'axe), avec $p=\cos(\theta/2)$ et
$||\overrightarrow{u}||=\sin(\theta/2)$.
$\theta\in[0,\pi]$, $\cos(\theta/2)$ est donc non négatif. Par
convention, $\sin(\theta)$ est non négatif. Le cas $\sin(\theta)$ négatif,
est équivalent au cas $\sin(\theta)$ positif en changeant $\overrightarrow{u}$ en $-\overrightarrow{u}$.

Le vecteur $\overrightarrow{u}$ n'est donc pas unitaire dans ces
formules. $q$ peut être manipulé comme un vecteur à 4 composantes
$(p,a,b,c)$ tout comme un complexe peut être manipulé comme un
vecteur à 2 composantes.

Liens avec la matrice de rotation $\mathcal{R}$ (voir formules d'Euler-Rodrigues) :
\[
\mathcal{R}=
\begin{pmatrix}
p^2+a^2-b^2-c^2 & 2(ab-cp) & 2(ac+bp) \\
2(ab+cp) & p^2-a^2+b^2-c^2 & 2(bc-ap )\\
2(ac-bp) & 2(bc+ap) & p^2-a^2-b^2+c^2\\
\end{pmatrix}
\]
La matrice $\mathcal{R}$ est une matrice orthogonale directe, ce qui signifie que ses colonnes forment une base orthonormée directe et donc son déterminant vaut 1.

Le produit de $\mathcal{R}$ par un vecteur $\overrightarrow{v}\left[\begin{array}[m]{l}
x\\
y\\
z
\end{array}\right]$ donne :

\[
\mathcal{R}\left[\begin{array}[m]{l}
x\\
y\\
z
\end{array}
\right] =
(2p^2-1)\left[\begin{array}[m]{l}
x\\
y\\
z
\end{array}\right]
+2p\left[\begin{array}[m]{l}
a\\
b\\
c
\end{array}\right]
\wedge
\left[\begin{array}[m]{l}
x\\
y\\
z
\end{array}\right]
+2(ax+by+cz)\left[\begin{array}[m]{l}
a\\
b\\
c
\end{array}\right]
\]
Avantages des quaternions sur les matrices:
\begin{enumerate}
  \item Coût calcul (et espace mémoire) moins élevé dans les produits de quaternions
(songez aux jeux videos ou on en enchaîne beaucoup).
  \item Pas de perte d'orthogonalité dans les produits de quaternions, il
suffit de renormaliser à 1 en cas de nécessité. Quand on
multiplie de nombreuses matrices de rotations, les arrondis
nous conduisent à une matrice qui ne vérifie plus ${}^t\mathcal{R}\mathcal{R}=I$
(matrice identité), et ce n'est pas facile à renormaliser.
\end{enumerate}

Avantage sur les angles d'Euler: c'est bien plus simple, plus
les avantages précédents !
\end{document}


 Voir aussi ::
http://fr.wikipedia.org/wiki/Quaternions_et_rotation_dans_l%27espace

Aucun commentaire:

Enregistrer un commentaire