Recherche d'un quaternion représentant la rotation d'un vecteur à un autre

105

J'ai deux vecteurs u et v. Existe-t-il un moyen de trouver un quaternion représentant la rotation de u vers v?

sdfqwerqaz1
la source

Réponses:

115
Quaternion q;
vector a = crossproduct(v1, v2);
q.xyz = a;
q.w = sqrt((v1.Length ^ 2) * (v2.Length ^ 2)) + dotproduct(v1, v2);

N'oubliez pas de normaliser q.

Richard a raison de dire qu'il n'y a pas de rotation unique, mais ce qui précède devrait donner «l'arc le plus court», ce dont vous avez probablement besoin.

Polaris878
la source
30
Sachez que cela ne gère pas le cas des vecteurs parallèles (à la fois dans la même direction ou pointant dans des directions opposées). crossproductne sera pas valide dans ces cas, vous devez donc d'abord vérifier dot(v1, v2) > 0.999999et dot(v1, v2) < -0.999999, respectivement, et soit renvoyer une identité quat pour les vecteurs parallèles, soit renvoyer une rotation de 180 degrés (autour de n'importe quel axe) pour les vecteurs opposés.
sinisterchipmunk
11
Une bonne implémentation de ceci peut être trouvée dans le code source d'ogre3d
João Portela
4
@sinisterchipmunk En fait, si v1 = v2, crossproduct serait (0,0,0) et w serait positif, ce qui se normalise à l'identité. Selon gamedev.net/topic/ ... cela devrait fonctionner correctement aussi pour v1 = -v2 et dans leur voisinage proche.
jpa
3
Comment quelqu'un a-t-il fait fonctionner cette technique? D'une part, sqrt((v1.Length ^ 2) * (v2.Length ^ 2))simplifie à v1.Length * v2.Length. Je ne pouvais obtenir aucune variation de ceci pour produire des résultats raisonnables.
Joseph Thomson
2
Oui, cela fonctionne. Voir le code source . L61 gère si les vecteurs font face à des directions opposées (retourne PI, sinon il renvoie l'identité par la remarque de @ jpa). L67 gère les vecteurs parallèles: mathématiquement inutiles, mais plus rapides. L72 est la réponse de Polaris878, en supposant que les deux vecteurs sont une longueur unitaire (évite un sqrt). Voir aussi les tests unitaires .
sinisterchipmunk
63

Solution de vecteur à demi-chemin

J'ai trouvé la solution que je crois qu'Imbrondir essayait de présenter (bien qu'avec une petite erreur, ce qui explique probablement pourquoi sinisterchipmunk avait du mal à la vérifier).

Étant donné que nous pouvons construire un quaternion représentant une rotation autour d'un axe comme ceci:

q.w == cos(angle / 2)
q.x == sin(angle / 2) * axis.x
q.y == sin(angle / 2) * axis.y
q.z == sin(angle / 2) * axis.z

Et que le point et le produit croisé de deux vecteurs normalisés sont:

dot     == cos(theta)
cross.x == sin(theta) * perpendicular.x
cross.y == sin(theta) * perpendicular.y
cross.z == sin(theta) * perpendicular.z

Étant donné qu'une rotation de u vers v peut être obtenue en tournant par thêta (l'angle entre les vecteurs) autour du vecteur perpendiculaire, il semble que nous puissions construire directement un quaternion représentant une telle rotation à partir des résultats du point et des produits croisés ; cependant, tel quel, thêta = angle / 2 , ce qui signifie que cela entraînerait le double de la rotation souhaitée.

Une solution consiste à calculer un vecteur à mi-chemin entre u et v , et à utiliser le point et le produit croisé de u et le vecteur à mi-chemin pour construire un quaternion représentant une rotation de deux fois l'angle entre u et le vecteur à mi-chemin , ce qui nous amène à v !

Il existe un cas particulier, où u == -v et un vecteur unique à mi-chemin devient impossible à calculer. Ceci est attendu, étant donné les rotations infiniment nombreuses de "l'arc le plus court" qui peuvent nous faire passer de u à v , et nous devons simplement tourner de 180 degrés autour de tout vecteur orthogonal à u (ou v ) comme solution de cas particulier. Ceci est fait en prenant le produit croisé normalisé de u avec tout autre vecteur non parallèle à u .

Le pseudo-code suit (évidemment, en réalité, le cas particulier devrait tenir compte des inexactitudes en virgule flottante - probablement en comparant les produits scalaires à un certain seuil plutôt qu'à une valeur absolue).

Notez également qu'il n'y a pas de cas particulier où u == v (le quaternion d'identité est produit - vérifiez et voyez par vous-même).

// N.B. the arguments are _not_ axis and angle, but rather the
// raw scalar-vector components.
Quaternion(float w, Vector3 xyz);

Quaternion get_rotation_between(Vector3 u, Vector3 v)
{
  // It is important that the inputs are of equal length when
  // calculating the half-way vector.
  u = normalized(u);
  v = normalized(v);

  // Unfortunately, we have to check for when u == -v, as u + v
  // in this case will be (0, 0, 0), which cannot be normalized.
  if (u == -v)
  {
    // 180 degree rotation around any orthogonal vector
    return Quaternion(0, normalized(orthogonal(u)));
  }

  Vector3 half = normalized(u + v);
  return Quaternion(dot(u, half), cross(u, half));
}

La orthogonalfonction renvoie tout vecteur orthogonal au vecteur donné. Cette implémentation utilise le produit croisé avec le vecteur de base le plus orthogonal.

Vector3 orthogonal(Vector3 v)
{
    float x = abs(v.x);
    float y = abs(v.y);
    float z = abs(v.z);

    Vector3 other = x < y ? (x < z ? X_AXIS : Z_AXIS) : (y < z ? Y_AXIS : Z_AXIS);
    return cross(v, other);
}

Solution Quaternion à mi-chemin

C'est en fait la solution présentée dans la réponse acceptée, et elle semble être légèrement plus rapide que la solution vectorielle à mi-chemin (~ 20% plus rapide par mes mesures, mais ne me croyez pas sur parole). Je l'ajoute ici au cas où d'autres comme moi seraient intéressés par une explication.

Essentiellement, au lieu de calculer un quaternion à l'aide d'un vecteur à mi-chemin, vous pouvez calculer le quaternion qui entraîne deux fois la rotation requise (comme détaillé dans l'autre solution), et trouver le quaternion à mi-chemin entre cela et zéro degré.

Comme je l'ai expliqué précédemment, le quaternion pour le double de la rotation requise est:

q.w   == dot(u, v)
q.xyz == cross(u, v)

Et le quaternion pour une rotation nulle est:

q.w   == 1
q.xyz == (0, 0, 0)

Le calcul du quaternion à mi-chemin est simplement une question de sommation des quaternions et de normalisation du résultat, tout comme avec les vecteurs. Cependant, comme c'est aussi le cas avec les vecteurs, les quaternions doivent avoir la même grandeur, sinon le résultat sera biaisé vers le quaternion avec la plus grande grandeur.

Un quaternion construit à partir du point et produit croisé de deux vecteurs aura le même ordre de grandeur que ces produits: length(u) * length(v). Plutôt que de diviser les quatre composants par ce facteur, nous pouvons à la place augmenter le quaternion d'identité. Et si vous vous demandez pourquoi la réponse acceptée complique apparemment les choses en utilisant sqrt(length(u) ^ 2 * length(v) ^ 2), c'est parce que la longueur au carré d'un vecteur est plus rapide à calculer que la longueur, nous pouvons donc enregistrer un sqrtcalcul. Le résultat est:

q.w   = dot(u, v) + sqrt(length_2(u) * length_2(v))
q.xyz = cross(u, v)

Et puis normalisez le résultat. Le pseudo code suit:

Quaternion get_rotation_between(Vector3 u, Vector3 v)
{
  float k_cos_theta = dot(u, v);
  float k = sqrt(length_2(u) * length_2(v));

  if (k_cos_theta / k == -1)
  {
    // 180 degree rotation around any orthogonal vector
    return Quaternion(0, normalized(orthogonal(u)));
  }

  return normalized(Quaternion(k_cos_theta + k, cross(u, v)));
}
Joseph Thomson
la source
12
+1: Génial! Cela a fonctionné comme un charme. Devrait être la réponse acceptée.
Rekin
1
La syntaxe du quaternion est activée sur certains exemples (Quaternion (xyz, w) et Quaternion (w, xyz)). Il semble également que dans le dernier bloc de code, les radians et les degrés soient mélangés pour exprimer les angles (180 vs k_cos_theta + k).
Guillermo Blasco
1
Quaternion (float, Vector3) est une construction à partir d'un vecteur scalaire, tandis que Quaternion (Vector3, float) est une construction à partir d'un axe-angle. Peut-être potentiellement déroutant, mais je pense que c'est correct. Corrigez-moi si vous pensez toujours que c'est faux!
Joseph Thomson
Ça a marché! Merci! Cependant, j'ai trouvé un autre lien similaire et bien expliqué pour effectuer l'opération ci-dessus. Je pensais que je devrais partager pour mémoire;)
pécheur
1
@JosephThomson La solution du quaternion à mi-chemin semble venir d'ici .
legends2k
6

Le problème tel qu'énoncé n'est pas bien défini: il n'y a pas de rotation unique pour une paire de vecteurs donnée. Prenons le cas, par exemple, où u = <1, 0, 0> et v = <0, 1, 0> . Une rotation de u à v serait une rotation pi / 2 autour de l'axe z. Une autre rotation de u vers v serait une rotation pi autour du vecteur <1, 1, 0> .

Richard Dunlap
la source
1
En fait, n'y a-t-il pas un nombre infini de réponses possibles? Parce qu'après avoir aligné le vecteur «de» avec le vecteur «vers», vous pouvez toujours faire tourner librement le résultat autour de son axe? Savez-vous quelles informations supplémentaires peuvent généralement être utilisées pour contraindre ce choix et rendre le problème bien défini?
Doug McClean
5

Pourquoi ne pas représenter le vecteur à l'aide de quaternions purs? C'est mieux si vous les normalisez d'abord peut-être.
q 1 = (0 u x u y u z ) '
q 2 = (0 v x v y v z )'
q 1 q rot = q 2
Pré-multipliez avec q 1 -1
q rot = q 1 -1 q 2
où q 1 -1 = q 1 conj / q norme
Cela peut être considéré comme une «division de gauche». La division droite, ce qui n'est pas ce que vous voulez, c'est:
q rot, right = q 2 -1 q 1

madratman
la source
2
Je suis perdu, la rotation de q1 à q2 n'est-elle pas calculée comme suit: q_2 = q_rot q_1 q_rot ^ -1?
yota
4

Je ne suis pas très bon sur Quaternion. Cependant, j'ai lutté pendant des heures à ce sujet et je n'ai pas pu faire fonctionner la solution Polaris878. J'ai essayé de pré-normaliser v1 et v2. Normaliser q. Normaliser q.xyz. Pourtant, je ne comprends toujours pas. Le résultat ne m'a toujours pas donné le bon résultat.

En fin de compte, j'ai trouvé une solution qui l'a fait. Si cela aide quelqu'un d'autre, voici mon code de travail (python):

def diffVectors(v1, v2):
    """ Get rotation Quaternion between 2 vectors """
    v1.normalize(), v2.normalize()
    v = v1+v2
    v.normalize()
    angle = v.dot(v2)
    axis = v.cross(v2)
    return Quaternion( angle, *axis )

Un cas particulier doit être fait si v1 et v2 sont parallèles comme v1 == v2 ou v1 == -v2 (avec une certaine tolérance), où je pense que les solutions devraient être Quaternion (1, 0,0,0) (pas de rotation) ou Quaternion (0, * v1) (rotation de 180 degrés)

Imbrondir
la source
J'ai une implémentation fonctionnelle, mais celle-ci est plus jolie, donc je voulais vraiment qu'elle fonctionne. Malheureusement, il a échoué tous mes cas de test. Mes tests ressemblent tous à quelque chose comme ça quat = diffVectors(v1, v2); assert quat * v1 == v2.
sinisterchipmunk
Il est peu probable que cela fonctionne du tout car il angletire sa valeur d'un produit scalaire.
sam hocevar
Où est la fonction Quaternion ()?
Juin Wang
3

Certaines des réponses ne semblent pas envisager la possibilité que le produit croisé soit égal à 0. L'extrait ci-dessous utilise la représentation de l'axe des angles:

//v1, v2 are assumed to be normalized
Vector3 axis = v1.cross(v2);
if (axis == Vector3::Zero())
    axis = up();
else
    axis = axis.normalized();

return toQuaternion(axis, ang);

Le toQuaternionpeut être implémenté comme suit:

static Quaternion toQuaternion(const Vector3& axis, float angle)
{
    auto s = std::sin(angle / 2);
    auto u = axis.normalized();
    return Quaternion(std::cos(angle / 2), u.x() * s, u.y() * s, u.z() * s);
}

Si vous utilisez la bibliothèque Eigen, vous pouvez également simplement faire:

Quaternion::FromTwoVectors(from, to)
Shital Shah
la source
toQuaternion(axis, ang)-> vous avez oublié de préciser ce que c'estang
Maksym Ganenko
Le deuxième paramètre fait anglepartie de la représentation axe-angle du quaternion, mesuré en radians.
Shital Shah
On vous a demandé de faire pivoter quaternion d'un vecteur à un autre. Vous n'avez pas d'angle, vous devez d'abord le calculer. Votre réponse doit contenir le calcul de l'angle. À votre santé!
Maksym Ganenko
C'est du C ++? qu'est-ce que ux ()?
Juin Wang
Oui, c'est C ++. u est le type de vecteur de la bibliothèque Eigen (si vous en utilisez un).
Shital Shah
2

Du point de vue de l'algorithme, la solution la plus rapide se présente sous forme de pseudocode

 Quaternion shortest_arc(const vector3& v1, const vector3& v2 ) 
 {
     // input vectors NOT unit
     Quaternion q( cross(v1, v2), dot(v1, v2) );
     // reducing to half angle
     q.w += q.magnitude(); // 4 multiplication instead of 6 and more numerical stable

     // handling close to 180 degree case
     //... code skipped 

        return q.normalized(); // normalize if you need UNIT quaternion
 }

Assurez-vous que vous avez besoin de quaternions unitaires (généralement, il est nécessaire pour l'interpolation).

REMARQUE: les quaternions non unitaires peuvent être utilisés avec certaines opérations plus rapidement que l'unité.

minorlogic
la source