Matrices de rotations intermédiaires

[neoirto]

Bonjour,

Je cherche un moyen de calculer 10 Matrices de rotations intermédiaires (Mrn) entre deux matrices de rotation 3x3 Ma et Mb, sous contraintes.
Il me faudra donc répondre aux différentes questions suivantes :
- Avec Ma et Mb connues : comment calculer la matrice moyenne (Mr) des 2 ?
- Avec Ma et Mb connues, comment calculer 10 matrices intermédiaires afin de décomposer le mouvement entre les deux matrices : par exemple : Mr1 = 0.1*Ma + 0.9*Mb, puis Mr2 = 0.2*Ma + 0.8*Mb... etc
- Comment ajouter une contrainte, comme par exemple fixer à une valeur prédéterminée un des angles d'euler sur toutes les matrices intermédiaires ?

Merci d'avance pour vos réponses !

[fatal_error]

bjr,

tes questions n'ont pas trop de sens à mon avis

comment calculer la matrice moyenne (Mr) des 2 ?

la moyenne (a,b) pour moi elle est symétrique, idem moyenne(a,b)==moyenne(b,a)

mais le produit matriciel est pas à priori commutatif A*B!=B*A sauf si A B sont particulières. Genre des matrices de rotation sur le même axe par ex..
Prends la somme des matrices... c'est faisable, mais je vois pas à quoi ca peut te servir.

Avec Ma et Mb connues, comment calculer 10 matrices intermédiaires afin de décomposer le mouvement

sachant que tu veux décomposer de Ma à Mb en 10 matrices intermédiraires, je m'attends à voir
Mb=M9*M8...*M0*Ma
or Mb contient n^2 coeffs, donc on a n^2 equations mais n^2*10 inconnues...

[Dlzlogic]

Bonjour,
Je vais vous décrire la méthode que l'utiliserais si j'avais à faire cela.
Mais je suppose que ce n'est pas la meilleure méthode.
En considérant qu'il n'y a pas de translation, la formule de rotation en 3D s'écrit ainsi :
X = XX.x + XY.y + XZ.z
Y = YX.x + YY.y + YZ.z
Z = ZX.x + ZY.y + ZZ.z
x, y, z sont les coordonnées de départ, X, Y, Z les coordonnées d'arrivée
XX, XY, ... ZZ sont les paramètres de la transformation.
Ces valeurs sont directement fonction des angles d'Euler, mais leur formulation est difficile, donc je ne les utilise pas..
Je pense que dans une première approche (expérimentale) pour passer de la position 1 à la position 2, il suffirait de faire passer les 9 paramètres, chacun indépendamment, de la position 1 à la position 2 par un pas proportionnel, c'est à dire 1/10 .
L'essai est très rapide à faire, cependant je ne suis pas sûr que ce soit bon.
La vérification consiste à calculer pour chaque position la distance au centre de rotation. Si celle-ci n'est pas bonne une simple homothétie devrait résoudre le problème.
Pour la contrainte relative à un angle, là, ça me parait plus compliqué. Je pense que j'essaierais de faire une grille de référence.

[leon1789]

@dlzlogic

Je pense que dans une première approche (expérimentale) pour passer de la position 1 à la position 2, il suffirait de faire passer les 9 paramètres, chacun indépendamment, de la position 1 à la position 2 par un pas proportionnel, c'est à dire 1/10 .

Si je comprends tes 10 lignes d'explication, tu proposes de calculer
Mr1 = 0.1*Ma + 0.9*Mb, puis Mr2 = 0.2*Ma + 0.8*Mb... etc
Mais cela ne fonctionne pas (sinon neoirto ne serait pas venu poser sa question) car les Mr ne sont pas des matrices de rotation.

@neoirto
vous vous intéressez à deux matrices de rotation dans R^3 (chacune caractérisée par un axe "orienté" de rotation et un angle "signé") : pour passer de l'une à l'autre, il me paraît naturel de passer progressivement l'axe de rotation de la première sur celui de la seconde, tout en modifiant progressivement l'angle de rotation de la première vers celui de la seconde.

Est-ce que cela paraît utilisable dans votre contexte ?

[Dlzlogic]

@dlzlogic

Si je comprends tes 10 lignes d'explication, tu proposes de calculer
Mr1 = 0.1*Ma + 0.9*Mb, puis Mr2 = 0.2*Ma + 0.8*Mb... etc
Mais cela ne fonctionne pas (sinon neoirto ne serait pas venu poser sa question) car les Mr ne sont pas des matrices de rotation.

C'est effectivement ce qu'il y a à craindre, ce sont des transformations affines (rotation + affinité), d'où une correction par homothétie.
De toute façon, il y a une infinité de solutions pour passer de Ma à Mb, il faut en trouver une, pas trop bête et rapide à appliquer, c'est à dire en évitant le calcul de lignes trigonométriques.

Sinon, il a toujours la solution, la première qui vient à l'esprit : calculer les paramètres en fonction des 3 angles et faire varier les 3 angles, à mon avis, c'est la solution lourde.

[neoirto]

Merci à tous pour vos réponses rapides !

vous vous intéressez à deux matrices de rotation dans R^3 (chacune caractérisée par un axe "orienté" de rotation et un angle "signé") : pour passer de l'une à l'autre, il me paraît naturel de passer progressivement l'axe de rotation de la première sur celui de la seconde, tout en modifiant progressivement l'angle de rotation de la première vers celui de la seconde.

@leon1789: absolument

@Dlzlogic:

il faut en trouver une, pas trop bête et rapide à appliquer, c'est à dire en évitant le calcul de lignes trigonométriques.

absolument !

Sinon, il a toujours la solution, la première qui vient à l'esprit : calculer les paramètres en fonction des 3 angles et faire varier les 3 angles, à mon avis, c'est la solution lourde.

Et bien, j'ai deja essayé :
- extraction des angles d'euler de Ma et Mb
- passage des pondérations aux angles d'euler
- vérifiage d'intervalle -PI // PI
- reconstruction de la matrice intermédiaire à partir des angles d'euler retravaillés comme expliqué ci-dessus
= Gros bricolage fonctionnant pas trop mal (probablement un peu lourd en calcul vu le nombre de fonctions trigo utilisées), mais surtout : à l'approche du GIMBAL LOCK : ça devient faux et ça tourne un peu dans tous les sens.

@fatal_error:
Oui j'ai bien conscience qu'à additionner des choux par des chèvres pondérées, je ne risque pas d'obtenir le graal ! Ma formulation initiale n'était là que pour "éclaircir" ma question. Ce n'est ni une base de réflexion mathématique, ni même comme du pseudo-code :
Mr1 = 0.1*Ma + 0.9*Mb, puis Mr2 = 0.2*Ma + 0.8*Mb = égale TOUT faux je suis bien d'accord !!!

Mb=M9*M8...*M0*Ma

Oui, c'est bien ce que je cherche à faire.
Ou encore reformulé : je cherche la matrice de passage (Mp) de Ma à Mb... Divisée par 10 (Mp10), et je ferai alors :
Mr1 = Mp10 * Ma;
Mr2 = Mp10 * Mr1;
...

or Mb contient n^2 coeffs, donc on a n^2 equations mais n^2*10 inconnues...

Oups, et bien j'espère qu'il existe une resolution trigo plus simple !

Sinon, concernant la contrainte : je pense qu'il suffirai de modifier Ma et Mb, de façon à ce qu'elles répondent toutes deux à la contrainte, et ensuite de calculer les matrice intermédiaires à partir de Ma' et Mb', non ?

Merci à tous !

[Dlzlogic]

Petite information supplémentaire.
La formulation de la transformation affine, translation + rotation + homothétie + affinité comporte 12 paramètres.
Ici, on veut faire abstraction de la translation, de l'homothétie et de l'affinité, et ne conserver que la rotation. A mon avis, ce n'est possible que sur le papier. Dans le cas présent, on travaille avec des nombres (coordonnées et angles connus par leurs valeurs trigo) et on fait des multiplications et des additions, il restera toujours un facteur correspondant à l'affinité qui ne sera pas 1 et à la translation qui ne sera pas 0.
Pour définir une rotation en 3D il faut théoriquement 3 points (non situés dans l'un des 3 plans de référence). De nombreux essais numériques ont montré que ça marchait pas. Les écarts sont faibles, mais non négligeables. La solution consiste à créer un 4è point pour former un tétraèdre à peu près bien équilibré, et on a ainsi 4 points qui permettent de calculer les 12 paramètres nécessaires.

Un autre point qu'il ne faut pas oublier : chacun des 9 paramètre peuvent passer par un maximum. Ci c'est le cas il faut calculer le maximum et modifier le pas de façon à atteindre ce maximum, puis s'en éloigner.

[neoirto]

@Dlzlogic:
Je pense ne pas avoir saisi l'histoire du maximum... Est-ce une histoire de minimisation du genre moindre carrés ?
Mais je confirme que je ne m'intéresse qu'à la rotation, donc hors homothétie et translation.

Petite question subsidiaire (toujours en langage "naturel" car je ne suis guère matheux) : comme le suggère leon1789, une matrice de rotation, ça peut aussi s'écrire comme un axe de rotation, et un angle de rotation autour de cet axe.
Je cherche donc la matrice de passage Mp telle que Mb = Mp * Ma

Ensuite j'extrait l'axe de rotation de la matrice, et je divise tout bêtement par 10 l'angle de rotation, avant de reconstruire Mp10 à partir de ce nouvel angle de rotation ???

Normalement, ça marche ça, non ?

Seulement, le coup des transferts Mat->axis-angle puis axis-angle->Mat, ça fait un paquet de racine à calculer... Existerait-il un moyen de simplifier les calculs puisque dans mon cas, je ne cherche à modifier que l'angle, et pas l'axe ?

[Dlzlogic]

Concernant le langage mathématique et en particulier l'usage du terme "matrice", son utilisation etc. je laisse cela à fatal-error et à Léon.
Dans le contexte des premiers messages, il semble bien qu"il s'agit de rotation en 3D de centre O. Si c'est une rotation par rapport à un axe, c'est beaucoup moins difficile, puisque c'est une rotation dans un plan autour d'un point.
Pour l'histoire du maximum, prenez le début de la formule
X=XX.x ...
XX a forcément un maximum = 1.
Il faut savoir, si on veut un joli rendu, que entre le début et la fin, on pourra avoir la composante X = x, donc XX est égal à 1 (les autres vivent leur vie de leur côté) dans l'intervalle début-fin. Mais si on fait varier XX proportionnellement entre le début et la fin, on ne passera forcément pas par son maximum XX=1.

Si vous êtes en 2D ou si on peut s'y ramener, alors, on peut utiliser les angles.
En fait, la rotation dans le plan départ -> arrivée est l'une des plus logiques, parce que la plus courte et mathématiquement, de loin la plus simple.
Donc, il s'agit de faire un changement de repère pour se trouver dans le plan voulu et faire une simple rotation en 2D.
Ce problème m'amuse, n'hésitez pas.

[leon1789]

Dans le contexte des premiers messages, il semble bien qu"il s'agit de rotation en 3D de centre O. Si c'est une rotation par rapport à un axe, c'est beaucoup moins difficile, puisque c'est une rotation dans un plan autour d'un point.

Depuis le début, neoirto parle de matrice 3x3, donc on est évidemment en 3D. Et en 3D, toute rotation possède un axe de rotation... Tu fais de l'humour ?

Encore une fois, c'est inutilement laborieux de penser avec 9 ou 12 paramètres.

[leon1789]

Ici, on veut faire abstraction de la translation, de l'homothétie et de l'affinité, et ne conserver que la rotation. A mon avis, ce n'est possible que sur le papier. Dans le cas présent, on travaille avec des nombres (coordonnées et angles connus par leurs valeurs trigo) et on fait des multiplications et des additions, il restera toujours un facteur correspondant à l'affinité qui ne sera pas 1 et à la translation qui ne sera pas 0.
Pour définir une rotation en 3D il faut théoriquement 3 points (non situés dans l'un des 3 plans de référence). De nombreux essais numériques ont montré que ça marchait pas. Les écarts sont faibles, mais non négligeables. La solution consiste à créer un 4è point pour former un tétraèdre à peu près bien équilibré, et on a ainsi 4 points qui permettent de calculer les 12 paramètres nécessaires.

Arrête donc tes histoires... C'est ta vie, ton expérience dans certains contextes. ICI, le contexte est clairement présenté et trouver une transformation convenable est a priori possible.

[leon1789]

Petite question subsidiaire (toujours en langage "naturel" car je ne suis guère matheux) : comme le suggère leon1789, une matrice de rotation, ça peut aussi s'écrire comme un axe de rotation, et un angle de rotation autour de cet axe.
Je cherche donc la matrice de passage Mp telle que Mb = Mp * Ma

Ensuite j'extrait l'axe de rotation de la matrice, et je divise tout bêtement par 10 l'angle de rotation, avant de reconstruire Mp10 à partir de ce nouvel angle de rotation ???

Normalement, ça marche ça, non ?

effectivement, je passerai de l'angle de Ma à celui de Mb par 10 pas constant, par exemple.
Disons angle_Mr1 = 0.1*angle_Ma + 0.9*angle_Mb

Seulement, le coup des transferts Mat->axis-angle puis axis-angle->Mat, ça fait un paquet de racine à calculer... Existerait-il un moyen de simplifier les calculs puisque dans mon cas, je ne cherche à modifier que l'angle, et pas l'axe ?

Vous vous voulez dire que les axes de vos deux rotations Ma et Mb sont les mêmes ? (ce qui ne peut que simplifier le problème !)

[neoirto]

Paix les amis ! Paix !

Vu qu'en fait c'est de ma faute et que vous avez probablement tous les deux raisons (du moins d'après ce que je comprends), je ne décris pas correctement mon problème depuis le début.

Je confirme que depuis le début, il s'agit bien de rotation 3D de centre O.

Et dans ma dernière réponse, je suggère de passer en 2D en utilisant une representation axis angle de LA MATRICE DE PASSAGE Mp, ce qui simplifierai les calculs tout en apportant la paix dans les foyers.
Une grosse astuce trigo qui évitera de balancer du Levenberg Marquardht (béni soit ce dernier pour l'orthographe de son patronyme). Mais est-ce bien possible ?

Sinon du coup je vois nettement mieux l'histoire des maximums, merci Dlzlogic.

[neoirto]

Vous vous voulez dire que les axes de vos deux rotations Ma et Mb sont les mêmes ? (ce qui ne peut que simplifier le problème !)

Béh non pour le coup !

avec la matrice de passage Mp telle que Mb = Mp * Ma, on en déduit Mp10eme construite à partir du meme axe que Mp et 1/10 * angle de Mp.

ensuite je ferais :
Mr1 = Mp10eme * Ma
Mr2 = Mp10eme * Mr1
Mr3 = Mp10eme * Mr2
Mr4 = Mp10eme * Mr3
Mr5 = Mp10eme * Mr4
Mr6 = Mp10eme * Mr5
Mr7 = Mp10eme * Mr6
Mr8 = Mp10eme * Mr7
Mr9 = Mp10eme * Mr8

et donc du coup :

Mb == Mp10eme * Mr9

Mais il y a peut être plus simple ??? Voire même c'est peut être faux comme raisonnement ?

[fatal_error]

re,

Si Mb et Ma sont des matrices de rotation, alors
Si on suppose Mp tq
Mb=MpMa
on peut ecrire
MbMa'=Mp car Ma est inversible et son inverse est sa transposée.
De fait, il reste à prouver que Mp est aussi une matrice de rotation.

Une matrice est dite de rotation si elle est
-orthogonale
-son determinant vaut 1

ici, det(Mp)=det(Mb)det(Ma')=1
Reste à montrer qu'elle est orthogonale.
De mémoire un groupe avec une loi c'est genre si on prend deux elements du groupe, et qu'on leur applique la loi alors le nouvel élément est toujours dans le groupe.
Si c'est bien le cas, comme Ma et Mb font parties du groupes des matrices orthogonales, alors Mp aussi et donc Mp est bien une matrice de rotation.

Du coup, elle fait juste un certain angle dans une certaine direction, et on peut prendre Mr=1/10Mp comme tu l'as dit.

En espérant ne pas écrire de boulettes

Edit: non j'ai encore débiter de l'error

Code: Tout sélectionner

%http://www.maths-forum.com/showthread.php?p=983361#post983361
a=pi/3;
Ma=[
1 0 0
0 cos(a) -sin(a)
0 sin(a) cos(a)
];
b=pi/6;
Mb=[
cos(b) 0 -sin(b)
0 1 0
sin(b) 0 cos(b)
];
%Mb=P*Ma
P=Mb*Ma';
r=10;
Mr=P/r;
U=Ma;
for i=1:r
  U=Mr*U;
end
shouldBeZeo=U-Mb %but is not! only for r==1




Code: Tout sélectionner


ans =

  -8.6603e-01   3.5083e-11   5.0000e-01
   3.5083e-11  -1.0000e+00  -8.9574e-11
  -5.0000e-01   8.9574e-11  -8.6603e-01


probablement le Mp/r qui est incorrect, on devrait s'en sortir avec la racine nieme de Mp, j'essaie de retrouver les trucs spectraux et compagnie..
avec

Code: Tout sélectionner

Mr=P^(1/r);
U=Ma;
for i=1:r
  U=Mr*U;
end
Mb-U


On a bien la matrice nulle

Code: Tout sélectionner

ans =

  7.7716e-16 - 1.6050e-16i  2.3580e-16 + 2.4830e-16i  -1.0547e-15 - 7.6502e-17i
  7.5322e-16 + 1.8479e-16i  -4.4409e-16 - 1.6457e-16i  4.1164e-16 + 2.1426e-16i
  -5.5511e-17 + 3.1669e-16i  1.3965e-15 - 5.2028e-16i  4.4409e-16 + 2.6568e-16i


il reste à montrer qu'on peut toujours prendre la racine nieme de Mp et que celle-ci est toujours une matrice de rotation.

[neoirto]

Cooool

Si c'est bien boulette-free, ce sera pas mal efficace, car avec les axis-angle, pas de de situation de GIMBAL LOCK, c'est bien ça ?

Bon ça fait beaucoup de racines à l'extraction et pas mal de multiplications en sens inverse, mais pourquoi pas si il n'y a pas plus simple.

Concernant la contrainte : par exemple pour fixer la rotation sur l'axe x du repere à zéro, je modifie Ma et Mb en conséquence et ça devrait marcher ?
Ca se fait comment cette modification svp, je veux dire au plus simple ?

EDIT: mince alors ça marche pas ??? En tout cas c'est sympa de m'aider avec tout ça, parce que je n'ai pas Matlab, je dois tout coder en C et je n'ai pas encore les transformations d'axis angle (il y a des singularités dans les transfos : voir liens ci dessus, ça vient peut être de là ?)

RE-EDIT: Ok, donc exit les racines et tous les calculs de mes liens ci-dessus : ça semble encore plus simple que prévu... D'enfer ton truc spectral là, même si je ne vois pas bien pourquoi la puissance.

[neoirto]

il reste à montrer qu'on peut toujours prendre la racine nieme de Mp et que celle-ci est toujours une matrice de rotation.

Aie, et c'est grave ?

Tout ce que je sais, c'est qu'on peut toujours prendre la transposée d'une matrice de rotation et calculer la matrice de passage comme tu l'as fait, "car je l'ai fait un paquet de fois sans jamais trouver de bug" (notez la qualité de démonstration).

Mais euu... Disons que le coté spectral de la puissance nième me passe un peu au-dessus là...

[leon1789]

il reste à montrer qu'on peut toujours prendre la racine nieme de Mp et que celle-ci est toujours une matrice de rotation.

Si j'ai bien compris ton idée, on pose Mp = Mb * inverse(Ma)
Il faut alors construire une matrice de rotation M telle que M^10 = Mp
et dans ce cas, on pose Mri = M^i * Ma .

Bon, il est clair que Mp est une matrice de rotation (car produit de deux matrices de rotations)
et donc il est facile de trouver une matrice M convenable : il suffit de trouver l'axe de rotation de Mp , puis l'angle de rotation .

Alors on pourrait prendre pour M la rotation d'axe le même que Mp et d'angle
Du coup, M^10 = Mp et chaque Mri = M^i * Ma est assurément une rotation.

[Ben314]

à l'approche du GIMBAL LOCK : ça devient faux et ça tourne un peu dans tous les sens.

Question c... : c'est quoi le GIMBAL LOCK ?
Sinon, comme tout le monde, je pese que le mieux, c'est de passer "continument" d'un axe à l'autre (en restant dans le plan qui contient les axes de départ et d'arrivé) et de passer "continument" d'un angle à l'autre.
Mais il y a déjà (donc avant de programmer) des petits problèmes techniques : le signe de l'angle d'une rotation de R^3 ne peut se définir qu'une fois l'axe orienté.
Par exemple, pour passer d'une rotation R1 d'axe orienté D et d'angle alpha à la rotation R2 de même axe D orienté dans l'autre sens et d'angle -alpha, il vaudrait mieux que ton programme ne fasse pas faire 180° à l'axe et change les angles : ça ferait con vu que R1=R2 !!!
Plus chiant, si R1 et R2 ont le même axe (orienté) et des angle opposés, on peut soit garder le même axe et modifier l'angle petit à petit, soit conserver l'angle et faire tourner petit à petit l'axe de 180°.
Lequel des deux préférer ? (et en fait tu aura systématiquement ce problème dans tout les cas vu qu'il y a deux axes orientés pour une même rotation...)

Dernière remarque : de façon théorique, quitte à multiplier (par exemple à gauche) toute tes matrices (y compris celles "intermédiaires" que tu cherche) par l'inverse de la rotation R_1, on peut parfaitement se contenter de traiter le cas où R_1=Id qui a l'énorme avantage... d'avoir absolument ce qu'on veut comme axe (par exemple le même que R_2...) Donc il y a juste à faire varier l'angle de 0 à celui de R_2 en gardant toujours le même axe.

Edit : ma "dernière remarque", ben c'est exactement ce que dit léon juste au dessus...

[leon1789]

Peut-être une idée pour réaliser facilement tous les calculs :

Mp = Mb * transpose(Ma) ;

t = l'angle de Mp (au signe près) = arccos( ( trace(Mp) - 1 ) / 2 )

Q = (Mp - transpose(Mp) ) / (2 sin(t))

P = (Mp/2 + transpose(Mp)/2 - cos(t)*I ) / (1-cos(t))

de sorte que Mp = cos(t) * I + sin(t) * Q + (1-cos(t)) * P

et la matrice M est donnée par

M = cos(t/10) * I + sin(t/10) * Q + (1-cos(t/10)) * P

M est bien une matrice de rotation telle que
Mri = M^i * Ma soit une matrice de rotation pour tout i
M^10 = Mp
Mb = M^10 * Ma

EDIT : pour des explications mathématiques, voir http://fr.wikipedia.org/wiki/Matrice_de_rotation#Axe_de_rotation

Ne pas oublier que
transpose(rotation) = rotation inverse = rotation d'angle opposé