Encontrando quatérnio representando a rotação de um vetor para outro

105

Eu tenho dois vetores uev. Existe uma maneira de encontrar um quatérnio que representa a rotação de u para v?

sdfqwerqaz1
fonte

Respostas:

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

Não se esqueça de normalizar q.

Richard está certo sobre não haver uma rotação única, mas o acima exposto deve fornecer o "arco mais curto", que é provavelmente o que você precisa.

Polaris878
fonte
30
Esteja ciente de que isso não lida com o caso de vetores paralelos (ambos na mesma direção ou apontando em direções opostas). crossproductnão será válido nesses casos, portanto, primeiro você precisa verificar dot(v1, v2) > 0.999999e dot(v1, v2) < -0.999999, respectivamente, e retornar um quat de identidade para vetores paralelos ou retornar uma rotação de 180 graus (em torno de qualquer eixo) para vetores opostos.
sinisterchipmunk
11
Uma boa implementação disso pode ser encontrada no código-fonte do ogre3d
João Portela
4
@sinisterchipmunk Na verdade, se v1 = v2, o produto cruzado seria (0,0,0) e w seria positivo, o que normaliza a identidade. De acordo com gamedev.net/topic/… deve funcionar muito bem também para v1 = -v2 e nas suas proximidades.
jpa
3
Como alguém fez essa técnica funcionar? Por um lado, sqrt((v1.Length ^ 2) * (v2.Length ^ 2))simplifica para v1.Length * v2.Length. Não consegui obter nenhuma variação disso para produzir resultados razoáveis.
Joseph Thomson
2
Sim, funciona. Veja o código fonte . L61 trata se os vetores estão voltados para direções opostas (retorne PI, caso contrário, retornaria identidade por observação de @jpa). L67 lida com vetores paralelos: matematicamente desnecessário, mas mais rápido. L72 é a resposta do Polaris878, assumindo que ambos os vetores têm comprimento unitário (evita um sqrt). Veja também testes de unidade .
sinisterchipmunk
63

Solução vetorial de meio caminho

Eu encontrei a solução que acredito que Imbrondir estava tentando apresentar (embora com um pequeno erro, provavelmente por isso que o sinistrochipmunk teve problemas para verificar).

Dado que podemos construir um quatérnio representando uma rotação em torno de um eixo assim:

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

E que o ponto e o produto cruzado de dois vetores normalizados são:

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

Visto que uma rotação de u para v pode ser alcançada girando por teta (o ângulo entre os vetores) em torno do vetor perpendicular, parece que podemos construir diretamente um quatérnion representando tal rotação a partir dos resultados dos produtos pontuais e cruzados ; no entanto, como está, theta = ângulo / 2 , o que significa que fazer isso resultaria no dobro da rotação desejada.

Uma solução é a de calcular um meio caminho vector entre u e v , e usar o ponto e o produto cruzado de u e a meio caminho vector para construir um Quatérnion representando uma rotação de duas vezes o ângulo entre u e a meio caminho vector, o que nos leva até v !

Há um caso especial, onde u == -v e um vetor intermediário único se torna impossível de calcular. Isso é esperado, dadas as infinitas rotações de "arco mais curto" que podem nos levar de u a v , e devemos simplesmente girar 180 graus em torno de qualquer vetor ortogonal a u (ou v ) como nossa solução de caso especial. Isso é feito tomando o produto vetorial normalizado de u com qualquer outro vetor não paralelo a u .

Segue-se um pseudocódigo (obviamente, na realidade, o caso especial teria que levar em conta imprecisões de ponto flutuante - provavelmente verificando os produtos escalares em relação a algum limite em vez de um valor absoluto).

Observe também que não há nenhum caso especial quando u == v (o quatérnio de identidade é produzido - verifique e veja por si mesmo).

// 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));
}

A orthogonalfunção retorna qualquer vetor ortogonal ao vetor fornecido. Esta implementação usa o produto vetorial com o vetor de base mais ortogonal.

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);
}

Solução de Quatérnio de Meio Caminho

Esta é realmente a solução apresentada na resposta aceita, e parece ser marginalmente mais rápida do que a solução vetorial intermediária (~ 20% mais rápida pelas minhas medições, embora não acredite apenas no que eu digo). Estou adicionando aqui para o caso de outras pessoas como eu estarem interessadas em uma explicação.

Essencialmente, em vez de calcular um quatérnio usando um vetor intermediário, você pode calcular o quatérnio que resulta em duas vezes a rotação necessária (conforme detalhado na outra solução) e encontrar o quatérnio a meio caminho entre aquele e zero graus.

Como expliquei antes, o quatérnio para o dobro da rotação necessária é:

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

E o quatérnio para rotação zero é:

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

Calcular o quatérnio a meio caminho é simplesmente uma questão de somar os quatérnios e normalizar o resultado, assim como acontece com os vetores. Porém, como também é o caso com vetores, os quatérnions devem ter a mesma magnitude, caso contrário, o resultado será desviado para o quatérnio com a magnitude maior.

A quaternion construído a partir do ponto e produto cruzado de dois vetores terá a mesma magnitude em que esses produtos: length(u) * length(v). Em vez de dividir todos os quatro componentes por esse fator, podemos aumentar o quaternion de identidade. E se você está se perguntando por que a resposta aceita aparentemente complica as coisas com o uso sqrt(length(u) ^ 2 * length(v) ^ 2), é porque o comprimento ao quadrado de um vetor é mais rápido de calcular do que o comprimento, portanto, podemos salvar um sqrtcálculo. O resultado é:

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

E então normalize o resultado. Segue pseudocódigo:

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
fonte
12
+1: Ótimo! Isso funcionou como um encanto. Deve ser a resposta aceita.
Rekin
1
A sintaxe do quaternion é ativada em alguns exemplos (Quaternion (xyz, w) e Quaternion (w, xyz)). Também parece que no último bloco de código, radianos e graus são misturados para expressar ângulos (180 vs. k_cos_theta + k).
Guillermo Blasco
1
Quaternion (float, Vector3) é a construção do vetor escalar, enquanto o Quaternion (Vector3, float) é construído a partir do ângulo do eixo. Talvez potencialmente confuso, mas acho que está correto. Corrija-me se ainda achar que está errado!
Joseph Thomson
Funcionou! Obrigado! No entanto, encontrei outro link semelhante e bem explicado para realizar a operação acima. Achei que deveria compartilhar para registro;)
pecador
1
@JosephThomson A solução do quatérnio a meio caminho parece vir daqui .
legends2k
6

O problema conforme afirmado não está bem definido: não há uma rotação única para um determinado par de vetores. Considere o caso, por exemplo, onde u = <1, 0, 0> e v = <0, 1, 0> . Uma rotação de u para v seria uma rotação pi / 2 em torno do eixo z. Outra rotação de u para v seria uma rotação pi em torno do vetor <1, 1, 0> .

Richard Dunlap
fonte
1
Na verdade, não existe um número infinito de respostas possíveis? Porque depois de alinhar o vetor "de" com o vetor "para", você ainda pode girar livremente o resultado em torno de seu eixo? Você sabe quais informações extras normalmente podem ser usadas para restringir essa escolha e tornar o problema bem definido?
Doug McClean
5

Por que não representar o vetor usando quatérnions puros? É melhor se você normalizá-los primeiro, talvez.
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é-multiplique com q 1 -1
q rot = q 1 -1 q 2
onde q 1 -1 = q 1 conj / q norma
Isso pode ser considerado "divisão à esquerda". A divisão certa, que não é o que você quer, é:
q rot, right = q 2 -1 q 1

madratman
fonte
2
Estou perdido, não é a rotação de q1 para q2 calculada como q_2 = q_rot q_1 q_rot ^ -1?
yota
4

Não sou muito bom no Quaternion. No entanto, eu lutei por horas nisso e não consegui fazer a solução Polaris878 funcionar. Eu tentei pré-normalizar v1 e v2. Normalizando q. Normalizando q.xyz. Mesmo assim, não entendo. O resultado ainda não me deu o resultado certo.

No final, porém, encontrei uma solução que o fez. Se isso ajudar mais alguém, aqui está meu código de trabalho (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 )

Um caso especial deve ser feito se v1 e v2 são paralelos como v1 == v2 ou v1 == -v2 (com alguma tolerância), onde acredito que as soluções devem ser Quaternion (1, 0,0,0) (sem rotação) ou Quaternion (0, * v1) (rotação de 180 graus)

Imbrondir
fonte
Eu tenho uma implementação funcionando, mas esta sua é mais bonita, então eu realmente queria que funcionasse. Infelizmente, falhou em todos os meus casos de teste. Todos os meus testes são parecidos quat = diffVectors(v1, v2); assert quat * v1 == v2.
sinisterchipmunk
É improvável que isso funcione, já que angleobtém seu valor de um produto escalar.
sam hocevar
Onde está a função Quaternion ()?
junho Wang
3

Algumas das respostas não parecem considerar a possibilidade de que o produto cruzado possa ser 0. O snippet abaixo usa a representação do eixo do ângulo:

//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);

O toQuaternionpode ser implementado da seguinte forma:

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);
}

Se você estiver usando a biblioteca Eigen, também pode simplesmente fazer:

Quaternion::FromTwoVectors(from, to)
Shital Shah
fonte
toQuaternion(axis, ang)-> você se esqueceu de especificar o que éang
Maksym Ganenko
O segundo parâmetro é o angleque faz parte da representação do ângulo do eixo do quaternion, medido em radianos.
Shital Shah
Você foi solicitado a fazer o quaternion girar de um vetor para outro. Você não tem ângulo, você tem que calculá-lo primeiro. Sua resposta deve conter o cálculo do ângulo. Felicidades!
Maksym Ganenko
Isso é c ++? o que é ux ()?
junho Wang
Sim, isso é C ++. u é o tipo de vetor da biblioteca Eigen (se você estiver usando uma).
Shital Shah
2

Do ponto de vista do algoritmo, a solução mais rápida olha em pseudocódigo

 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
 }

Certifique-se de que você precisa de quatérnios de unidade (normalmente, é necessário para interpolação).

NOTA: Quatérnions não unitários podem ser usados ​​com algumas operações mais rápidas do que a unidade.

minorlogic
fonte