かみのメモ

コンピュータビジョン・プログラムな話題中心の勉強メモ

回転ベクトル・回転行列・クォータニオン・オイラー角についてまとめてみた

以前の記事(OpenCVで取得したカメラパラメータをUnityで使う - かみのメモ)を書いたのをきっかけに、三次元座標系での回転の表現方法について色々調べたのでまとめておきたいと思います。

はじめに

三次元座標系で回転を表現するための方法として、回転ベクトル, 回転行列, オイラー, クォータニオン四元数がよく知られています。

この記事ではこれら4つの表現方法について

  1. 原理とその特徴
  2. 右手系・左手系の変換
  3. 各表現の相互変換(代表的なもののみ)

の3つを紹介していきます。

実際にPythonで回転後の座標を計算したり各表現を相互変換したりするプログラムは以下の記事で紹介しています↓。

なおこの記事はコンピュータビジョンと航空力学をかじっただけの人が書いたものです。 できるだけ誤りのないように書いているつもりですが、もし間違いを見つけた場合は報告していただけるとありがたいです。

スポンサーリンク

もくじ

0. 準備

0.1. 回転を表現するとは?

はじめに回転を表現するとはどういうことなのかを復習しておきます。

工学の分野では三次元空間にあるオブジェクトの姿勢を考えるとき、どこか(例えば地面)に基準となる世界座標系を、オブジェクトに固有のローカル座標系を設定してその関係を回転成分並進成分で記述するという考え方をします。

この回転成分の記述方法としては、以下4つの表現のどれかを使うことがほとんどです。

  1. 世界座標系で見ると (x,y,z) である点は、ローカル座標系で見ると (x,y,z) になる
  2. ローカル座標系で見ると (u,v,w) である点は、世界座標系で見ると (u,v,w) になる
  3. 世界座標系の座標軸をこう回転させると、ローカル座標系の座標軸に重なる
  4. ローカル座標系の座標軸をこう回転させると、世界座標系の座標軸に重なる

これらの表現は表裏一体で、どれで記述してもオブジェクトの回転すなわち姿勢を表現していることになります。

今回紹介する回転ベクトル回転行列クォータニオンは1,2の発想に近いものです。 "変換前の座標系での座標値"から"変換後の座標系での座標値"への射影を関数で記述する、つまり点を別の点に移動させるツールで回転を定義しようというアイデアです。

一方、オイラーは3,4の発想に近いもので、"変換前の座標系の座標軸"と"変換後の座標系の座標軸"の関係をなるべく直感的に記述しようというアイデアに基づいています。

回転を考えるときは自分が1〜4のどれをベースに考えているのかということを常に意識しましょう。 特に1と3、2と4は似て非なる発想です(むしろ1と4、2と3の方が近い)。 この辺りを曖昧にしたまま直感的に考えると失敗しやすいので注意が必要です。

0.2. 左手系と右手系の関係

次に右手座標系と左手座標系の関係について決めごとをしておきます。

右手座標系と左手座標系は軸のどれか1本を反転させた関係であると言えます。 べつに右手系X,Y,Zが左手系X,-Z,-Yと対応するような変換を考えてもいいのですが、右手系X,Y,Zと左手系X,-Y,Zの対応を考えた後にX軸周りに270°回転させたものと考えたほうが分かりやすいです。 ということで、この記事では右手系と左手系は2本の軸を共有していて1本の軸だけ反転している関係と考えます。

今後この記事で右手・左手変換の例を紹介するときはY軸反転、つまり右手系で見ると (1,2,3) である点を左手系で見ると (1,2,3) になるような変換を取り上げることにします。

f:id:kamino-dev:20190728202742p:plain:w600
Y軸反転の右手/左手座標系の関係

1. 回転ベクトル

まずは回転ベクトル(rotation vectorの紹介です。

回転ベクトルは「ある点 p を、原点を通るベクトル v の周りでその長さ |v| の分の角度だけ回転させる」という形で回転(=点の移動)を記述したものです。

f:id:kamino-dev:20190728203207p:plain:w270
回転ベクトルの概念図(右手座標系)

たとえば (π/2,0,0) はX軸を中心に90°回すように点を移動させることを意味します。 回転方向は、右手座標系では右ねじの方向(ベクトル方向を向いて時計回り)を正とし、左手座標系では左ねじの方向(ベクトル方向を向いて反時計回り)を正とします。

回転ベクトルはあくまで「点の回転移動」だけを記述するツールなので、グローバル->ローカルなのかローカル->グローバルなのかによって意味が変わることに注意してください。 つまり p, p のどちらがグローバル座標でどちらがローカル座標なのかは状況によって違うということです。 これは後で紹介する回転行列・クォータニオンも同じです。

1.1. 回転ベクトルの特徴

コンパクト

回転ベクトルは回転を最もコンパクトに記述する、つまり最低限の変数(3自由度)で記述する表現方法です。 そのためパラメータの自由度を拘束する条件式を考える必要がなく扱いやすいです。

ジンバルロックがない

後で紹介するオイラー角と違ってジンバルロックを生じないという特徴があります。 これは回転ベクトルを滑らかに変化させることで、ある回転の状態(=姿勢)から任意の別の回転の状態へと滑らかに変化させられることを意味します。

最適化問題で利用しやすい

上2つの特徴は、姿勢に関する最適化問題を解くときに非常に重要になります。

最適化問題数値計算で解くときは、各パラメータを少しずつ変化させて評価関数の値がよくなるところを探索するという方法を採ります。 このとき、パラメータに無駄な自由度があるとそれを抑えるための制約条件を加える必要が出てきますが、制約付き最適化問題は制約なしに比べて複雑になるので計算コストが高くなってしまいます。 また、パラメータを滑らかに変化させたときに実際の姿勢も滑らかに変化してくれないと、効率的に最適化を行うことができません。

このような理由から最適化問題では姿勢を回転ベクトル表現で記述することが多いです。 例えばOpenCVcalibrateCamera()solvePnP()なども内部では回転ベクトル表現を利用しています。

同じ回転を表す組

長さを ±2π だけ変化させた回転ベクトルは同じ回転を表していることになります。

回転量に対して 0|v|π という条件を課すことで、任意の回転は一意な回転ベクトルと一対一に対応します(ただし |v|=π のときだけ vx>0 など方向を一意に定める条件が必要)。

1.2. 回転ベクトルの逆変換

回転ベクトルは向きを逆にしてやれば、真逆の回転を表す回転ベクトルになります。 そのため (u,v,w) の逆変換は (u,v,w) になります。

1.3. ロドリゲスの回転公式

回転ベクトルを使って点を移動させるアルゴリズムロドリゲスの回転公式としてよく知られています。 回転ベクトル v によって点 p を回転させるとき、回転後の点 p は以下の式で表されます。

v¯=v / |v|θ=|v|p=pcosθ+v¯(v¯p)(1cosθ)+(v¯×p)sinθ

この式は右手系と左手系で共通です。 気になる方は v=(u,v,w), p=(x,y,z) を代入したとき(右手系)と、v=(u,v,w), p=(x,y,z) を代入したとき(Y軸反転の左手系)を比較して、p のY成分の符号だけが反転することを確かめてみてください。

1.4. 回転ベクトルの右手系と左手系の変換

右手系で取得した回転ベクトルを左手系での表現に直すときは、向きを反転させた後で普通のベクトルと同じように変換してやればよいです。 たとえば、右手系で観測した回転ベクトル (u,v,w) はY軸反転した左手系で観測すると (u,v,w) になります。

向きを反転させるのは、右手系と左手系で正とする回転方向が逆であるためです。

2. 回転行列

次にみんな大好き回転行列(rotation matrix)のお話です。

回転行列は三次元の回転を3×3の直交行列で表現する方法です。 点の座標を縦ベクトルで表現する場合、回転行列 R によって点 p を回転させると回転後の点 p は以下の式で表されます。

p=Rp

2.1. 回転行列の特徴

回転を9パラメータで表現する

回転行列の各列ベクトルは、変換前の座標系における基底ベクトルXYZが、変換後にどのベクトルに対応するかを示しています。 見た目上は3自由度の回転を9つのパラメータで表現することになりますが、直交行列であるという条件( RRT=RTR=E もしくは同値な条件として、列ベクトル表現を R=[r1 r2 r3] として r1×r2=r3,  r1r2=0,  |r1|=|r2|=1 )によって6自由度が抑えられるため、結果として3自由度を過不足なく表現できます。

あとで説明するように、回転行列はある回転ベクトルに対するロドリゲスの回転公式を成分に書き下したときの表現に相当します。

任意の回転と一対一に対応する

回転行列は任意の回転と一対一に対応します。

ジンバルロックがない

回転ベクトルと同じくジンバルロックは生じません。

数学的にシンプルに記述できる

回転行列表現の利点は、ベクトルの回転(点の移動)や回転の合成(重ねがけ)をベクトルと行列の式で記述できるようになることです。

たとえば座標系AからBへの座標値の変換とBからCへの座標値の変換がそれぞれ回転行列 BRA, CRB であるなら、座標系AからCへの座標値の変換は CRA=CRBBRA となりますし、座標系CからAへの座標値の変換は ARC=CRA1=CRAT=(CRBBRA)T となり、かなりシンプルに記述できます。

また回転行列を拡張した同次座標系の概念を導入すれば、ベクトルと行列の式で回転・並進移動をまとめて扱えるようになります。 幾何関係を数学的に分析したいときによく使われる表現方法です。

2.2. 回転行列の逆変換

上でちらっと取り上げたように、回転行列の逆行列は逆変換に相当します。 また直交行列の性質として R1=RT が成立するので、結局、転置行列 RT が逆変換を表す回転行列となります。

2.3. 回転行列と回転ベクトルの関係

先ほどの回転ベクトル vロドリゲスの回転公式を成分に書き下すと、

[pxpypz]=[cosθ+v¯x2(1cosθ)v¯xv¯y(1cosθ)v¯zsinθv¯xv¯z(1cosθ)+v¯ysinθv¯yv¯x(1cosθ)+v¯zsinθcosθ+v¯y2(1cosθ)v¯yv¯z(1cosθ)v¯xsinθv¯zv¯x(1cosθ)v¯ysinθv¯zv¯y(1cosθ)+v¯xsinθcosθ+v¯z2(1cosθ)][pxpypz]

となります。 ここで出てきた3×3の行列がそのまま回転行列になります。

ロドリゲスの回転公式は右手系・左手系で共通ですので、もちろんこの式も右手系・左手系で共通です。 ただし次の項で説明するように、右手系と左手系どちらで見るかによって v¯ の値が変わるので、実際の回転行列の成分値は異なります。

2.4. 回転行列の右手系と左手系の変換

回転行列の右手系から左手系への変換は上の式を調べればわかります。

例えば右手系で v=(u,v,w) である回転ベクトルは、Y軸反転の左手系では v=(u,v,w) になります。 これを上の式に代入すれば、Y軸反転の左手系での表現に変換するためには回転行列の1行2列目・2行1列目・2行3列目・3行2列目の4箇所の符号を反転させれば良いことがわかります。

3. クォータニオン

続いて最近流行り(?)のクォータニオン(quaternion : 四元数について紹介します。

3.1. クォータニオンの定義

たまに「クォータニオン=回転」みたいな表現を耳にしますが、クォータニオン自体は複素数を拡張したような数学的な体系の一つです。

名前の通り、クォータニオンq=w+xi+yj+zk のように4つの基底を持つベクトルとして表現されます。 そして乗算に関して次のようなルールを持っています。

i2=j2=k2=ijk=1

この式を分解して表にすると次のようになります。

× 1 i j k
1 1 i j k
i i 1 k j
j j k 1 i
k k j i 1

この表に従ってクォータニオン同士の積を成分に書き下すと以下のようになります。

q1q2=(w1w2x1x2y1y2z1z2)+(w1x2+x1w2+y1z2z1y2) i+(w1y2x1z2+y1w2+z1x2) j+(w1z2+x1y2y1x2+z1w2) k

試してみるとわかりますが、クォータニオンでは結合則が成立します( q1(q2q3)=(q1q2)q3 )。 一方、交換則は成立しません( q1q2q2q1 )。

また、あるクォータニオン q=w+xi+yj+zk に対して、q=wxiyjzk を共役クォータニオンと呼びます。

本当はもう少し色々と定義すべきことがあるのですが、回転を表すツールとして使う分にはここまでの知識で十分です。

3.2. クォータニオンによる回転の表現

クォータニオンによる回転表現は回転ベクトルとの関係を示すのが一番わかりやすいです。

回転ベクトル v=(a,b,c)(方向ベクトル v¯=(a¯,b¯,c¯) 、長さ θ=|v|)があるとき、同じ回転を表現するクォータニオン

q=cosθ2+ia¯sinθ2+jb¯sinθ2+kc¯sinθ2

であると定義されます。 この定義から w2+x2+y2+z2=1 の制約条件が付き、結果的に3自由度の回転をぴったり表現できていることになります。

上記の通りクォータニオンは4変数で表現されるので、本質的には4次元空間上での射影を考える必要があるのですが、3次元の回転を考えるだけならそこまでの分析は必要ないでしょう。 とりあえずは上の定義に従って回転ベクトルに読み替えるのが図解しやすく直観的にもわかりやすいと思います。

クォータニオンと回転行列の相互変換については、この定義と2.3. 回転ベクトルと回転行列の関係を照らし合わせてみてください。

クォータニオンを使って点 p=(x,y,z) を回転させるときは、点の座標を pq=0+xi+yj+zk というクォータニオンに読み替えて

pq=qpqq

を計算します。 pqi,j,k 成分を取り出すと、回転後の点の座標値になっています。

またクォータニオン q1 で回転させた後に q2 で回転させたような回転は q2q1 となります。 この回転にしたがって点を移動させるときは pq=(q2q1)pq(q2q1) のようにすればよいです。 ちなみに (q2q1)=q1q2 が成り立ちます。

3.3. クォータニオンの特徴

ジンバルロックがない

回転ベクトルと同じくジンバルロックは生じません

同じ回転を表す組

クォータニオン qq は同じ回転を表します。 これは、方向 v¯ 長さ θ の回転ベクトルと、方向 v¯ 長さ 2πθ の回転ベクトルが同じ回転を表すことに対応しています。 そのため最終的な回転の状態としては同じですが、右回りで辿り着いたのか、左回りで辿り着いたのか、という区別があります。

クォータニオンの成分 (w,x,y,z) を順に見ていったとき最初に出てくる非ゼロな成分が正である、という条件を付ければクォータニオンと任意の回転は一対一に対応します。

球面線形補間が楽に計算できる

回転をクォータニオンで表現することのメリットとして球面線形補間(slerp : spherical linear interpolation)が楽に計算できることが挙げられます。 球面線形補間はある回転の状態(=姿勢)から別の回転の状態への変化を角度に対して線形に補間する処理です。 つまり、オブジェクトをある姿勢から別の姿勢までスムーズに回転させるような処理が簡単に実装できるようになります。

回転を表す2つのクォータニオン q1, q2 と配分率 0h1 が与えられたとき、補間して得られるクォータニオンは以下のように計算されます。

qslerp=q1(q1q2)h=q2(q2q1)1h=(q2q1)hq1=(q1q2)1hq2

累乗の計算はコストが高いので以下のように近似することもあるようです。

cosΩ=q1q2q~slerp=q1sin((1h)Ω)+q2sinhΩsinΩ

またslerpの際に q2 の代わりに q2 を使うと補間されるルートが変わります。 ざっくり言うと、一方は地球を西廻りして日本→韓国のルートを辿りますが、もう一方は日本→ブラジル→アフリカ→韓国という遠回りの東廻りルートを辿ります。 短い方のルートを選びたいときは、q1内積を取って正になる方を使えばよいです。

3.4. クォータニオンの逆変換

クォータニオン q の逆変換は、その共役クォータニオン q となります。

3.5. クォータニオンの右手系と左手系の変換

クォータニオンの右手系・左手系の変換は、反転していない2軸に対応する成分の符号を反転してやればよいです。 気になる方はクォータニオンによる回転の定義と回転ベクトルの右手系・左手系の変換を照らし合わせてみてください。

例えば右手系で q=w+xi+yj+zk であったクォータニオンは、Y軸反転の左手系では q=wxi+yjzk となります。

4. オイラー

続いてオイラー角(Euler angles)の紹介です。 これが一番ややこしいです。

オイラー角は「座標軸の1本を選びその周りで回転させる、という操作を3回繰り返す」という形で座標系の関係を表現する方法です。 3回それぞれの回転量を並べて (θ1,θ2,θ3) の形で回転を表現します。

オイラー角と一口に言っても、選ぶ回転軸の組み合わせによってX->Y->Z, Z->Y->X, Z->Y->Z, X->Y->Xなど12通りの表現があります。 ただ、結局はどの姿勢を (0,0,0) と置くかの違いなので、どの表現も基本的な性質は変わりません。

また2回目・3回目の回転のときに「もとの座標系の軸」と「回転後の座標系の軸」のどちらを基準に回転させるのかによっても種類が分かれます。 回転前の世界座標系の軸を使うものを外因性オイラー角(extrinsic Euler angles)、回転後の座標系の軸を使うものを内因性オイラー角(intrinsic Euler angles)と呼ぶようです。 なおこの記事では便宜上、内因性・外因性と呼ぶことにしますが正式な日本語訳があるのかはよくわかりません。

一般的にオイラー角といえば内因性オイラー角を指すことがほとんどです。 また、4.2. 回転行列との関係で説明しますが、内因性オイラー角と外因性オイラー角には対応関係があり、基本的な性質は共通しています。

ちなみにX->Y->Zオイラー角をz-y-xとかx-y'-z"のように記述している文献もたまに見かけるので、文献ごとの流儀をよく確認しましょう。

4.1. オイラー角の特徴

リンク機構と対応している

内因性オイラー角を説明するときのイメージとしてはジンバルがよく取り上げられます。 これはジャイロスコープなどに採用されている機構で、外側のリングから順番に回転量を指定していくと中のコマの姿勢が1つに定まる、という機構がオイラー角と対応しています。 また、ロボットアームのようなリンク機構もオイラー角と対応しています。

f:id:kamino-dev:20191210084542p:plain:w300
ジンバル(From : https://upload.wikimedia.org/wikipedia/commons/e/e2/3D_Gyroscope.png

このように実在する回転機構と一致すること、回転をイメージしやすいことから、オイラー角は広く利用されています。 Unityのシーンエディタでもオブジェクトの回転を入力するのにX->Y->Zの内因性オイラー角が使用されています。

ヨー・ピッチ・ロールと対応している

航空を始めとするいくつかの工学分野では姿勢を「ヨー(yaw)、ピッチ(pitch)、ロール(roll)」の3パラメータで表現することがありますが、これもオイラー角と対応しています。

ただしこの表現はオイラー角と同じように、XYZを前後左右上下のどれに割り当てるかと、外因性オイラー角と内因性オイラー角のどちらを使うかによって微妙に定義が変わります。

よく使われるのは下をZ・右をY・前をXとしたZ->Y->Xの内因性オイラー角における回転量をヨー・ピッチ・ロールと呼ぶスタイルらしいですが、文献によってバラバラの可能性があるのでよく確認してください。

ジンバルロックがある

オイラー角の厄介な性質としてジンバルロック(特異姿勢)が知られています。 詳しい解説は他文献に譲りますが、これはある姿勢において回転の自由度が2に制限されてしまう現象です。 同じ意味ですが、3変数のうちの2つが同じ1自由度の表現に使われるため姿勢を一意に記述できなくなる現象、と説明されることもあります。

ジンバルロック時は、オイラー角を滑らかに変化させても、ある姿勢から別の姿勢へと滑らかに変化させられない場合があります。 これは画像から姿勢を推定するなどの最適化問題を解くときや、CGでオブジェクトを自由に回転させるときに問題になることがあります。

同じ回転を表現する組

回転角を 0θ2π の範囲で変化させるとき、同一の角度を表現するオイラー角は複数存在します。 特にジンバルロックの姿勢では1方向の回転を2軸が担っている状態なので、無数のオイラー角表現が存在します。 詳細に知りたい場合は、対象とするオイラー角を回転行列に変換して各成分を調べてみてください。

4.2. オイラー角の逆変換

オイラー角の逆変換は簡単で、回転の操作を向きを逆にして逆順に施してやればよいです。

例えばX->Y->Zオイラー(θx,θy,θz) の逆変換は、Z->Y->Xオイラー(θz,θy,θx) になります。 これをX->Y->Zオイラー角で表現するのはなかなか面倒なので省略します。

4.3. オイラー角と回転行列の関係

ここでは比較的よく使うオイラー角から回転行列への変換だけを紹介します。

まず座標系Aから座標系Bへの回転を表すX->Y->Zの内因性オイラー角の場合。 各軸周りの回転行列を Rx, Ry, Rz とすると全体の回転行列は ARB=RzRyRx となります。 座標系AからBへのオイラー角を考えていたのに、出てくるのが「座標系Bの座標値をAで見たときの座標値に変換する回転行列」である辺りが割と直感に反しているので気をつけてください。

回転量を (θx,θy,θz) として成分に書き下すと以下のようになります。

ARB=[cosθzsinθz0sinθzcosθz0001][cosθy0sinθy010sinθy0cosθy][1000cosθxsinθx0sinθxcosθx]

この式は右手系・左手系で共通です。 X->Y->Xなどの他の内因性オイラー角も同じ要領で計算できます。

次に座標系Aから座標系Bへの回転を表すX->Y->Zの外因性オイラー角の場合。 これは計算式を細かく追うと面倒なので省略しますが、座標系Bの基底が座標系Aで見るとどのようになっているのかを考えることで ARB=RxRyRz であることが確認できます。 よく見ると、内因性オイラー角と行列を掛け合わせる順を逆にしただけの式であることがわかります。

つまり、座標系Aから座標系Bへの回転を表すX->Y->Zの内因性オイラー(θx,θy,θz) と、座標系Aから座標系Bへの回転を表すZ->Y->Xの外因性オイラー(θz,θy,θx) は全く同じ回転を表現しています。

4.4. オイラー角の右手系と左手系の変換

オイラー角の右手系・左手系の変換は、反転している軸以外の回転量の符号を逆にしてやればよいです。

例えば右手系でのX->Y->Zオイラー(θx,θy,θz) は、Y軸反転左手系ではX->Y->Zオイラー(θx,θy,θz) となります。

5. あとがきと参考文献

過去最長の記事になってしまいました。 図を作る体力が残らなかったのでとりあえずこのまま投稿します…。

誤植・追記希望等ありましたらコメントください。

この記事を作成するに当たって以下の文献を参考にさせていただきました。 ありがとうございました。


もしよければ↓の☆を1クリックお願いします!