クォータニオン(四元数)な話をするよん
はじめに
あ、「はじめに」は読まなくていいです(ただの私事)
最近、マイクラから離れてたのですがちょうど3DCGのレンダラを作ってるときに「回転処理の実装でクォータニオンを使用する」的なこと言ってたら、フレンドが何やらマイクラのitem_displayなどのエンティティのNBTにクォータニオンがあると言ってきたんですね。
そんなこと言われちゃぁ気になりますよ。
調べた。Wiki見たらアフィン変換とか書いてあるじゃないですか。こんなのガチで数学とかやってるコマンダーしか理解できない奴じゃないですか。
面白いしそれなりに困ってる人いたので解説しようと思ったわけです。
仕入れたばかりの新しくて新鮮なクォータニオンの説明です
あ、あと割とマイクラに則した内容となります
必要になる数学的知識(たぶん)
- 三角関数
- ベクトル
- 簡単な数学(算数)
簡単な三角関数の説明
今回使うのはsinとcosと逆関数だけですのでsin,cosは斜辺の長さが1で一つの鋭角に対する直角三角形の比を計算する物だと思ってもらって構いません。
また、逆関数はsin,cos,tanの結果を角度に直すものだと思ってもらって構いません。
それぞれarcsin,arccos,arctanと表記します。
また、arctanに関しては原点から指定された座標に向けるような式を書くことができます。
ただ、atan2は例外で入力に座標を使います。本来arctanfracyxで行う入力ですが、これではある条件下で逆回転を示す値になってしまうのでそれを解決する関数です。
textatan2left(y,xright)と表記します。
簡単なベクトルの説明
NBTのMotionタグみたいな感じです。
n個の要素を持ち、それを同時に扱い、考えることができるやつです。
座標の概念をより数学的に使いやすくした物だと思って構わないです。
さっそくクォータニオン(四元数)ってなんやねん
さあクォータニオンって何なんでしょうねぇ？
数学では複素数の拡張として説明されますが今回はそんなこと理解する必要ないです。
クォータニオンとはズバリ...
3次元の全ての回転を4つの数で表した数！
数を数で表すとか言う最高に意味の分からない文章ですが、つまり4次元ベクトルで3次元のすべての回転を表しているとでも思ってください。
というかそっちの方がわかりやすいですね（笑）
ついでに、ここで言っておきますが、今後一つ一つの回転を行うクォータニオンを回転クォータニオン、そのエンティティと言った物体の傾きを表すクォータニオンを姿勢クォータニオンとしていうことがあります。
じゃあそのクォータニオンの数一つ一つは何を表しているのでしょうか？
クォータニオンのそれぞれの要素の意味
クォータニオンの要素のそれぞれの意味というのは実は回転と兼ねて考えることは非常に難しいものです。
なので今回は回転を表せる不思議な数として考えた方が手っ取り早いです
一応数学的にはそれぞれの要素を$x_0,x_1,x_2,x_3$とすると
x=x_0+x_1i+x_2j+x_3k\space \left(x_0,x_1,x_2,x_3\in\mathbb{R}\right)
となっています。
回転にクォータニオンを使用するメリットとデメリット
一番大きいメリットとしてジンバルロックが発生しないということが挙げられやすいです。
ですがそれ以外にも
- 実際に回転処理に組み込む際の計算量が圧倒的に少ない
- 軸指定ができる
- 球面回転補間が容易である
ということがあげられます。
逆にデメリットは何でしょうか？
一番多く取り上げられるのはその数値を見ただけではどれがどのような回転の要素になっているかイメージしにくいということです。
ほかにはないですかね
まあ唯一あるデメリットは仕組みを理解しにくいでしょうね
回転の要素とクォータニオンの要素の関係
回転の軸を$\left(a_x,a_y,a_z\right)$の単位ベクトルとし、その周りを時計回りに$\theta$ラジアン回転したとするとそれに対して次のように$\left(q_x,q_y,q_z,q_w\right)$($q_w+q_xi+q_yj+q_zk$)というクォータニオンに変換できます。
$q_x = a_x\sin\frac{\theta}{2}$
$q_y = a_y\sin\frac{\theta}{2}$
$q_z = a_z\sin\frac{\theta}{2}$
$q_w = \cos\frac{\theta}{2}$
割とシンプルですね。
また、クォータニオンと軸をx-y-zの順で$\alpha$、$\beta$、$\gamma$ラジアンと変化させるオイラー角についての変換は次の通りです。1
(逆はこちらの記事をご覧ください)
\alpha = \text{atan2}\left(2\left(q_xq_y+q_zq_w\right),q^2_x-q^2_y-q^2_z+q^2_w\right)
\beta = \arcsin\left(2\left(q_xq_z-q_yq_w\right)\right)
\gamma = \text{atan2}\left(2\left(q_xq_w+q_yq_z\right), q^2_x+q^2_y-q^2_z-q^2_w\right)
位置ベクトルとクォータニオンの計算
簡単な積の形で表せます
位置ベクトルを$\vec{v}$、回転クォータニオンを$q$、$q$とは逆回転(後述)の共役クォータニオンを$\bar{q}$とするとき、回転後の新しい位置ベクトル$\vec{v'}$は次のように表せる2
\vec{v'} = q\vec{v}\bar{q}
また、回転前の姿勢クォータニオン(どのくらいどの方向に傾いているかを表すクォータニオン)を$p$とし、回転後の姿勢クォータニオン$p'$は次のように表せる
p'=qp
クォータニオンの計算と注意点
さあ立て続けに計算式が出てきていますがクォータニオンはどのように計算する物なのでしょうか？
ちなみに同じ軸となる要素同士で掛け算してはいけません！(先ほどの計算式が成立しません)
その計算方法を説明します
簡単です。
任意のクォータニオン$p=\left(p_x,p_y,p_z,p_w\right)$と回転クォータニオン$q=\left(q_x,q_y,q_z,q_w\right)$に対して、その積を次のように定義する。
qp=\left(q_wp_x-q_zp_y+q_yp_z+q_xp_w, q_zp_x+q_wp_y-q_xp_z+q_yp_w, -q_yp_x+q_xp_y+q_wp_z+q_zp_w, -q_xp_x-q_yp_y-q_zp_z+q_wp_w\right)
計算式はそれこそ長いですが、すべてコマンドでも可能なくらい簡単な計算ですね！
3次元位置ベクトルとクォータニオンの積
3次元ベクトルとクォータニオンの計算は、3次元ベクトルのそれぞれの成分に加えて値が0のw成分をクォータニオンとして計算する。
計算結果のw成分は無視して3次元ベクトルとして使用する。
位置ベクトルを$\vec{v}=\left(x,y,z\right)$とすると
q\vec{v}\rightarrow q\vec{v'}=q\left(x,y,z,0\right) \left(\vec{v'}\in\mathbb{H}\right)
注意点1 積の交換法則は成り立たない
クォータニオンにおいて、$qp$と$pq$は等しくなりません。
注意点2 回転クォータニオンのノルムは必ず1
回転クォータニオンのノルム、つまり次の計算式で得られる結果は必ず1でなければなりません。
|q|=q_x^2+q_y^2+q_z^2+q_w^2
注意点3 同じ回転を示すクォータニオンの存在
$q$と$-q=\left(-q_x, -q_y, -q_z, -q_w\right)$は同じ回転を示します。
また、その亜種的なものとして逆回転を示す共役クォータニオンが存在します。
$q$に対して$\bar{q}$と書き、計算式は
\bar{q}=\left(-q_x, -q_y, -q_z, q_w\right)
です。
おわりに
多分コマンドで使用する分にはこれくらいを理解していれば十分なのではないでしょうか？
そして今回のこの記事で何回か目に見かけたであろう$i,j,k$ですが、この3つの数には次のような性質があります。
i^2=j^2=k^2=ijk=-1
本来はこれがあって初めて四元数となるのですが、プログラミングではそのように定義せずとも結局演算方法を定義する羽目になるので4次元ベクトルとしてとらえてよかったというわけです。
さて、私には理解できていないexecute幾何学とやらがコマンドの世界にはあるようですがそれで三角関数を作ることができるのでしょうか？
そんな動画を見たことあるので作れるのでしょうけども
詳しく知りたい人は四元数やらクォータニオンやらで調べてください()
私は力尽きました
参考文献
クォータニオン (Quaternion) を総整理！ ～ 三次元物体の回転と姿勢を鮮やかに扱う ～ - Qiita
おまけ-実装例(JavaScript)
頭が痛くなるような実装例。
動作確認してないなんて言えない
常に立方体の頂点がくるくる回転するプログラムです。(一言も画面上でとは言っていない)
// @ts-check
class Vector3 {
type = "Vector3";
/**
* Define constructor.
* @param {Number} x
* @param {Number} y
* @param {Number} z
*/
constructor(x, y, z) {
if (typeof(x) !== "number") throw TypeError();
if (typeof(y) !== "number") throw TypeError();
if (typeof(z) !== "number") throw TypeError();
this.x = x;
this.y = y;
this.z = z;
}
/**
* Is this object Vector3 class?
* @param {*} Object
* @returns {boolean}
*/
static IsVector3(Object) {
try {
if (Object.type == "Vector3") return true;
else return false;
}
catch (e) {
return false;
}
}
/**
* Define Vector3 + Vector3.
* @param {Vector3} v0
* @param {Vector3} v1
* @returns {Vector3}
*/
static add(v0, v1) {
return new Vector3(v0.x + v1.x, v0.y + v1.y, v0.z + v1.z);
}
/**
* Return Vector3's length.
* @returns {Number}
*/
Length() {
return this.x**2 + this.y**2 + this.z**2;
}
/**
* Define Vector3 to Quaternion.
* @returns {Quaternion}
*/
ToQuaternion() {
return new Quaternion(this.x, this.y, this.z, 0);
}
}
class Quaternion {
type = "Quaternion";
/**
* Define constructor.
* @param {Number} x
* @param {Number} y
* @param {Number} z
* @param {Number} w
*/
constructor(x=0, y=0, z=0, w=1) {
if (typeof(x) !== "number") throw TypeError();
if (typeof(y) !== "number") throw TypeError();
if (typeof(z) !== "number") throw TypeError();
if (typeof(w) !== "number") throw TypeError();
this.x = x;
this.y = y;
this.z = z;
this.w = w;
}
/**
* Is this object Quaternion class?
* @param {*} Obj
* @returns {boolean}
*/
static IsQuaternion(Obj) {
try {
if (Obj.type == "Quaternion") return true;
else return false;
}
catch (e) {
return false;
}
}
/**
* Trasform arbitrary axis Angle to Quaternion.
* @param {Number} angle
* @param {Vector3} AxisVector
* @returns {Quaternion}
*/
AngleAxis(angle, AxisVector) {
if (typeof(angle) !== "number") throw TypeError();
if (!Vector3.IsVector3(AxisVector)) throw TypeError();
const UnitVecAxis = new Vector3(AxisVector.x / AxisVector.Length(), AxisVector.y / AxisVector.Length(), AxisVector.z / AxisVector.Length());
const RadAngle = Math.PI * angle / 180;
this.x = UnitVecAxis.x * Math.sin(RadAngle / 2);
this.y = UnitVecAxis.y * Math.sin(RadAngle / 2);
this.z = UnitVecAxis.z * Math.sin(RadAngle / 2);
this.w = Math.cos(RadAngle / 2);
return this;
}
/**
* Define Quaternion * Quaternion.
* @param {Quaternion} q0
* @param {Quaternion} q1
* @returns {Quaternion}
*/
static multiply(q0, q1) {
if (!Quaternion.IsQuaternion(q0)) throw TypeError();
if (!Quaternion.IsQuaternion(q1)) throw TypeError();
return new Quaternion(
q0.w * q1.x - q0.z * q1.y + q0.y * q1.z + q0.x * q1.w,
q0.z * q1.x + q0.w * q1.y - q0.x * q1.z + q0.y * q1.w,
-q0.y * q1.x + q0.x * q1.y + q0.w * q1.z + q0.z * q1.w,
-q0.x * q1.x - q0.y * q1.y - q0.z * q1.z + q0.w * q1.w
);
}
/**
* Return conjugated Quaternion.
* @returns {Quaternion}
*/
Inverse() {
return new Quaternion(-this.x, -this.y, -this.z, this.w);
}
/**
* Define Quaternion to Vector3
* @returns {Vector3}
*/
ToVector3() {
return new Vector3(this.x, this.y, this.z);
}
/**
* Define Norm with Quaternion.
* @returns {Number}
*/
Norm() {
return this.x**2 + this.y**2 + this.z**2 + this.w**2;
}
}
/**
* Define rotational processing.
* @param {Vector3} Pos
* @param {Quaternion} RotationQua
* @returns {Vector3}
*/
function CalculateRotation(Pos, RotationQua) {
if (!Vector3.IsVector3(Pos)) throw TypeError();
if (!Quaternion.IsQuaternion(RotationQua)) throw TypeError();
let Rotation = new Quaternion();
if (RotationQua.Norm() != 0) {
Rotation = new Quaternion(RotationQua.x / RotationQua.Norm(), RotationQua.y / RotationQua.Norm(), RotationQua.z / RotationQua.Norm(), RotationQua.w / RotationQua.Norm());
}
return Quaternion.multiply(Quaternion.multiply(Rotation, Pos.ToQuaternion()), Rotation.Inverse())
.ToVector3();
}
// e.g.
function RenderClear(Scbuff) {
if (!Array.isArray(Scbuff)) throw TypeError();
for (let i=0; i<480; i++) {
for (let j=0; j<640; j++) {
Scbuff[i][j] = 0x000000;
}
}
}
let Screenbuffer = [];
const obj = {
points: [
new Vector3(-1, -1, -1),
new Vector3( 1, -1, -1),
new Vector3( 1, 1, -1),
new Vector3(-1, 1, -1),
new Vector3(-1, -1, 1),
new Vector3( 1, -1, 1),
new Vector3( 1, 1, 1),
new Vector3(-1, 1, 1)
],
pos: new Vector3(0, 0, 0),
rotation: new Quaternion()
};
const FocalLen = 1 / (2 * Math.tan(Math.PI / 2));
while (1) {
RenderClear(Screenbuffer);
obj.rotation = Quaternion.multiply(new Quaternion().AngleAxis(1, new Vector3(0, 1, 0)), obj.rotation);
obj.points.forEach((elem)=>{
let Pointpos = CalculateRotation(elem, obj.rotation);
let x = Math.floor(FocalLen * Pointpos.x * 200 / Pointpos.z);
let y = Math.floor(FocalLen * Pointpos.y * 200 / Pointpos.z);
Screenbuffer[y][x] = 0xFFFFFF;
});
}
/*
Copyright (c) 2024 @Scr_MIYUKINNGU2 MIYUKINNGU
*/