OpenGL/Kvaternio

Suunnan ja pyörimisliikkeen (eli rotaation) esittämiseen on useita vaihtoehtoja. Yaw, pitch ja roll -esitysmuoto aiheuttaa gimbal lock -ongelman. Rotaatiomatriisin jatkuva kertominen aiheuttaa epätarkkuuksia, jolloin se joudutaan joko muodostamaan uudestaan tai normalisoimaan silloin tällöin; molemmat operaatiot ovat hitaita.

Kvaternio (engl. quaternion) on kätevä esitysmuoto pyörimisliikkeelle. Se koostuu neljästä luvusta. Kvaternioja voidaan ketjuttaa kertomalla kuten matriiseja, epätarkkuus on helppo korjata ja niitä voidaan interpoloida yksinkertaisesti. Matemaattisesti kvaternio on hyvin abstrakti rakenne, eikä siihen kannattane paneutua kuin ehkä geometrian tai fysiikan kannalta; niissä sitä kutsutaan myös Eulerin parametreiksi.

Mielenterveyden säilyttämiseksi voi olla parasta käyttää kvaternioja sellaisenaan, ”taikakaluna”, ainakin aluksi.

Käytäntöjä

Tässä tekstissä kvaternio ilmaistaan neliulotteisena vektorina:

q = (q_x, q_y, q_z, q_w)

q on kvaternio ja q_x, q_y, q_z ja q_w ovat reaalilukuja.

Tässä kirjassa kvaternioista puhuttaessa tarkoitetaan aina rotaatiota esittävää kvaterniota. Rotaatio on aina yksikkökvaternio eli sen pituus on 1. Tästä poikkeaminen aiheuttaa virheitä; liukulukujen epätarkkuudesta johtuvaa poikkeamaa hoidetaan normalisoinnilla.

C-kielen esitystapa

C-kielessä kvaternio esitetään yleensä 4 alkion taulukkona tai tietueena. Tässä kirjassa käytetään tietuetta selkeyden vuoksi.

/* olkoon w viimeinen luku (indeksi 3) */
typedef float Quat[4];

tai

typedef struct {
    float x,y,z,w;
} Quat;

Joskus vektori- ja skalaariosa halutaan erottaa:

typedef float Vec[3];
typedef struct {
    Vec vector;
    float scalar;
} Quat;

C++:ssa ja muissa oliokielissä tehdään tietysti kvaternioluokka. Esimerkkikoodia on esimerkiksi NeHe Productionsin materiaalissa.

Muodostaminen

q muodostetaan kiertokulmasta θ ja akselin suunnan määräävästä yksikkövektorista (e_x, e_y, e_z) seuraavasti:

q_x = e_x sin(θ / 2)

q_y = e_y sin(θ / 2)

q_z = e_z sin(θ / 2)

q_w = cos(θ / 2)

Huomataan, että kiertämätön kvaternio on (0, 0, 0, 1). Joskus sitä kutsutaan identiteettikvaternioksi.

#include <math.h>
/* kiertokulma annetaan radiaaneina toisin kuin glRotatessa */
void muodosta_kvaternio_rotaatiovektorista(Quat* kohde, float kiertokulma, float x, float y, float z)
{
    float puoli_kiertokulmaa = kiertokulma * 0.5f;
    float sin_puoli_kiertokulmaa = sin(puoli_kiertokulmaa);
    kohde->x = x*sin_puoli_kiertokulmaa;
    kohde->y = y*sin_puoli_kiertokulmaa;
    kohde->z = z*sin_puoli_kiertokulmaa;
    kohde->w = cos(puoli_kiertokulmaa);
}

Kvaternion muodostaminen yaw, pitch ja roll -arvoista onnistuu kertomalla kolme kvaterniota halutussa järjestyksessä. Järjestys on päinvastainen kuin glRotate-käskyissä.

q = q_yaw q_pitch q_roll

Kertominen

Kvaternioiden kertolaskulla voidaan ketjuttaa rotaatioita matriisien tapaan. Järjestys on täysin päinvastainen! Viimeiseksi kerrottu kvaternio siis ”suoritetaan” viimeiseksi.

Olkoon q = ab. Tällöin:

q_x = a_w b_x + a_x b_w + a_y b_z − a_z b_y

q_y = a_w b_y + a_y b_w + a_z b_x − a_x b_z

q_z = a_w b_z + a_z b_w + a_x b_y − a_y b_x

q_w = a_w b_w − a_x b_x − a_y b_y − a_z b_z

/* tulos != a ja tulos != b */
void kerro_kvaterniot(Quat* tulos, const Quat* a, const Quat* b)
{
    tulos->x = a->w*b->x + a->x*b->w + a->y*b->z - a->z*b->y;
    tulos->y = a->w*b->y + a->y*b->w + a->z*b->x - a->x*b->z;
    tulos->z = a->w*b->z + a->z*b->w + a->x*b->y - a->y*b->x;
    tulos->w = a->w*b->w - a->x*b->x - a->y*b->y - a->z*b->z;
}

Muuntaminen rotaatiomatriisiksi

Jotta pisteitä voisi pyöritellä tehokkaasti, kvaternio muutetaan matriisiksi. Sitten sillä kerrotaan OpenGL:n nykyinen matriisi.

Olkoon kvaternio (x, y, z, w). Sitä vastaava rotaatiomatriisi M määritellään seuraavasti:

$\mathbf {M} =\left({\begin{array}{cccc}1-(2y^{2}+2z^{2})&2xy+2zw&2xz-2yw&0\\2xy-2zw&1-(2x^{2}+2z^{2})&2yz+2xw&0\\2xz+2yw&2yz-2xw&1-(2x^{2}+2y^{2})&0\\0&0&0&1\end{array}}\right)$

/* matriisi tuotetaan OpenGL:n haluamaan muotoon (pystyrivit peräkkäin taulukkoon) */
void matriisi_kvaterniosta(float* M, const Quat* q)
{
    float xx = q->x * q->x;
    float xy = q->x * q->y;
    float xz = q->x * q->z;
    float xw = q->x * q->w;
    float yy = q->y * q->y;
    float yz = q->y * q->z;
    float yw = q->y * q->w;
    float zz = q->z * q->z;
    float zw = q->z * q->w;

    M[ 0] = 1.0f - 2.0f * (yy + zz);
    M[ 1] =        2.0f * (xy - zw);
    M[ 2] =        2.0f * (xz + yw);

    M[ 4] =        2.0f * (xy + zw);
    M[ 5] = 1.0f - 2.0f * (xx + zz);
    M[ 6] =        2.0f * (yz - xw);

    M[ 8] =        2.0f * (xz - yw);
    M[ 9] =        2.0f * (yz + xw);
    M[10] = 1.0f - 2.0f * (xx + yy);

    M[3] = M[7] = M[11] = M[12] = M[13] = M[14] = 0.0f;
    M[15] = 1.0f;
}
...
float matriisi[16];
matriisi_kvaterniosta(matriisi, &kvaternio);
glMultMatrixf(matriisi);

Operaatio on suhteellisen nopea, koska monimutkaisin laskutoimitus on kertolasku.

Normalisointi

Jatkuvia kertolaskuja kärsinyt kvaternio on parempi normalisoida, jotta se pysyisi numeerisesti vakaana. Normalisointi tapahtuu kuten tavallisessa vektorissa: jokainen komponentti jaetaan pituudella.

float kvaternion_pituus(const Quat* q)
{
    return sqrt(q->x*q->x + q->y*q->y + q->z*q->z + q->w*q->w);
}
void normalisoi_kvaternio(const Quat* q)
{
    float kerroin = 1.0f / kvaternion_pituus(q);
    q->x *= kerroin;
    q->y *= kerroin;
    q->z *= kerroin;
    q->w *= kerroin;
}

Huomautuksia

• Tässä kvaternioista kerrottiin vain ohjelmoijalle välttämätön. Englanninkielisessä Wikipediassa on aiheesta runsaasti tarkempaa tietoa.
• Yaw, pitch ja roll -arvoista voitaneen muodostaa kvaternio suoraan, mikä on hitusen nopeampaa.
• Pyörimisliikkeen säilyttäminen matriisina saattaa olla nopeampaa, jos matriisi normalisoidaan vain harvoin. Toisaalta monimutkainen, kiinteän kappaleen dynamiikkaa (engl. rigid body dynamics) ratkova fysiikanmallinnus voi joutua kiertämään kappaleita hyvinkin usein.
• Usein käytetty tapa normalisoida rotaatiomatriisi on Gramin–Schmidtin menetelmä.
• OGRE-grafiikkaengine ja ODE-fysiikkaengine käyttävät kvaterniosta esitysmuotoa, jossa taulukon parametrit ovat järjestyksessä w, x, y, z.