[Mathematics] 21. Quaternion Slerp (Spherical Linear Interpolation)

Lerp 또는 Linear Interpolation은 예전 포스트에서 다룬 적이 있다.

간단히 말해 어떤 값 A와 B 사이를 P * t로 표현하는 법으로, 위 식을 보면 t=0일때 A이고 1일때 B임을 알 수 있다.

(A와 B는 스칼라 값 또는 벡터값 모두가 될 수 있다)

 

Quaternion의 경우에도 이러한 interpolation이 가능하다. 즉 어떤 두 사원수 q와 r 사이를 p라는 사원수로 선형으로 왔다갔다 할 수 있다는 것이다. (구면을 따라 움직이는 회전 운동을 Interpolate하는 것이므로 Spherical interpolation이라고 표현한다.)

 

Quaternion interpolation 혹은 slerp 공식을 구하는 법은 lerp와 기본적인 느낌은 비슷하다.

우리가 구하려는 식은 t가 0일때 q이고 t가 1일때 r인, p에 대한 식을 찾는 것이다.

즉 정의역은 t이고 치역은 p이며 치역의 범위는 q~r이다.

 

우선, d * q = r이 되는 사원수 d가 존재함은 자명하다.
(회전운동q를 회전운동r로 변환할 수 있는 회전운동d가 존재함은 자명하기 때문이다. 참고: https://gamesmith.tistory.com/161) 

고로

d * q = r

d = r * q^-1

 

이제 p = d q 꼴의 식에 t를 정의역으로 사용해 p를 q에서 r 사이로 왔다갔다 하게 만들어보자.

단순히 생각하면 p = dt * q를 하면 될 것 같지만 사원수 d에 스칼라 t를 곱하면 d의 크기가 변하므로 이는 옳지 않다.

(우리는 사원수의 크기를 1로 유지하고 싶다)

 

그렇다면 어떻게 해야 할까?

사원수의 구조를 살펴보자.

n hat 은 normal vector이다. 즉 회전하려는 축 방향으로 크기가 1인 벡터이다.

사원수는 다음과 같은 형태인데, 어떤 두 회전을 interpolate하려면 회전각을 interpolate해야 한다는 느낌을 직관적으로 받을 수 있다. 즉 d의 세타들에 t를 곱해주는 것으로 사원수를 interpolate할 수 있다.

이때 "d의 세타들에 t를 곱해주는 것"은 놀랍게도 d^t와 값이 동일하다.

즉 d의 세타들에 t를 곱한 값과 d에 d를 t번 곱해준 값은 같다.

 

Jorge의 마지막 말을 토대로 보았을 때 사원수의 제곱이 저런 형태로 계산되는 이유는 수체계가 변하면서 연산을 계산하는 방법도 변했기 때문인 듯 하다.

이렇게 정의내렸을 때, 우리는 q와 r을 interpolate하는 식을 구할 수 있다.

우리가 원하던 대로 p값이 t=0일때 q이고, t=1일때 r임을 알 수 있다.

 

코드로 나타내면 다음과 같다.

// Raising this quaternion to a power of t gives us the
// fraction t of the rotation of the quaternion
const Quaternion Quaternion::operator^(float t) const
{
	// Convert the quaternion back into axis/angle
	float a;
	Vector n;

	ToAxisAngle(n, a);

	// Scale the angle by t
	float at = a*t;

	// Make a new quaternion out of the new value
	return Quaternion(n, at);
}

// Spherical linear interpolation of a quaternion
const Quaternion Quaternion::Slerp(const Quaternion& other, float t) const
{
	const Quaternion& q = *this;
	Quaternion r = other;

	// 본 포스트에서 구한 공식으로 계산해도 잘 작동한다
	// return ((r * q.Inverted()) ^ t) * q;

	// 그러나 더 빠른 방식이 존재하는데, 이 공식을 도출하는 방법은 "3D Primer For Graphics
	// and Game Development" by Dunn and Parberry의 section 10.4.13를 참고
	float flCosOmega = w*r.w + r.v.Dot(v);
	if (flCosOmega < 0)
	{
		// Avoid going the long way around.
		r.w = -r.w;
		r.v = -r.v;
		flCosOmega = -flCosOmega;
	}

	float k0, k1;
	if (flCosOmega > 0.9999f)
	{
		// Very close, use a linear interpolation.
		k0 = 1-t;
		k1 = t;
	}
	else
	{
		// Trig identity, sin^2 + cos^2 = 1
		float flSinOmega = sqrt(1 - flCosOmega*flCosOmega);

		// Compute the angle omega
		float flOmega = atan2(flSinOmega, flCosOmega);
		float flOneOverSinOmega = 1/flSinOmega;

		k0 = sin((1-t)*flOmega) * flOneOverSinOmega;
		k1 = sin(t*flOmega) * flOneOverSinOmega;
	}

	// Interpolate
	Quaternion result;
	result.w = q.w * k0 + r.w * k1;
	result.v = q.v * k0 + r.v * k1;

	return result;
}

 

참고: https://youtu.be/fRSaaLtYj68?list=PLW3Zl3wyJwWOpdhYedlD-yCB7WQoHf-My