Mathematics for Computer Graphics
1
ユークリッド空間の変換
1.1
ベクトル積
2 つのベクトル a = (ax, ay, az), b = (bx, by, bz)に対し、ベクトル c = (cx, cy, cz)が、a, b に直交し，ノ
ルムが a, b の作る平行四辺形の面積に一致し、a, b, c の順に正の方向が x, y, z 軸の順に正の方向 (右手系) であるとき、 c を a, b のベクトル積 (vector product) または外積 (outer product) といい，a× b と表す．
c = ï¯ ¯¯ ¯ ay az by bz ¯¯ ¯¯ ¯, ¯¯ ¯¯ ¯ az ax bz bx ¯¯ ¯¯ ¯, ¯¯ ¯¯ ¯ ax ay bx by ¯¯ ¯¯ ¯ ! このとき (a× b, c) = ¯¯ ¯¯ ¯¯ ¯ ax ay az bx by bz cx cy cz ¯¯ ¯¯ ¯¯ ¯ (1) cの向きは、(1) の最後の式の値が正になる方向でも定められる。
1.2
直交行列
長さを変えない線形変換を直交変換という。 定義と定理 (1) 行列 U が Ut U = Eのとき，直交行列 (orthogonal matrix) と呼ぶ． (2) 行列 U が直交行列である同値な条件は、任意のベクトル a, b に対し、(U a, U b) = (a, b) であるこ とである。 (3) U が直交行列である同値な条件は U−1=t U . (4) U = (u1, u2, u3)が直交行列である同値な条件は{u1, u2, u3} が正規直交系 ((ui, uj) = δij)をなすこ とである．(5) U が直交行列のとき det(U ) = ±1 である．det(U) = 1 のとき回転行列 (rotation matrix) であり， det(U ) =−1 のときは鏡映 (reflection) を含む． (6) 直交行列の標準形は T (θ)±= ±1 0 0 0 cos θ − sin θ 0 sin θ cos θ である．+ のときは det(U ) = 1 となり，固有値 1 に対する固有ベクトルが回転軸になる．− のとき は det(U ) =−1 となり，固有値 −1 に対する固有ベクトルが鏡映面の法線ベクトルになる．
(7) 回転軸が u1，回転角が θ の回転行列 A は，{u1, u2, u3} を正規直交系とするとき，U = (u1, u2, u3) とおくと A = U T (θ)+ t U すなわち t U AU = T (θ)+ (2) 座標軸の回りの回転 座標軸を回転軸にして θ 回転の行列は，次の通りである。左から x–軸， y–軸，z–軸の回転を表す。 X(θ) = 1 0 0 0 cos θ − sin θ 0 sin θ cos θ , Y (θ) = cos θ 0 sin θ 0 1 0 − sin θ 0 cos θ , Z(θ) = cos θ − sin θ 0 sin θ cos θ 0 0 0 1 任意の軸での回転行列 初めに正規直交基底を求めるシュミットの正規直交化定理から始める． シュミットの正規直交化 {a1, a2, . . . , an, . . .} を一次独立なベクトルの列とすると，{a1, a2, . . . , an} の張 る空間内で，次の過程で正規直交基底{u1, u2, . . . , un| (ui, uj) = δij} の列を求めることができる． u1= a1/ p (a1, a1), k= 2 に対して tk = ak− Pk−1 i=1(ai, ui)ui, uk = tk/ p (tk, tk) 証明 1. u1= a1/ p (a1, a1)と置けば, (u1, u1) = 1 2. t2= a2− tu1 が u1に直交するために 0 = (t2, u1) = (a2, u1)− t(u1, u1) = (a2, u1)− t したがって t = (a2, u1). u2= t2/ p (t2, t2)と置けば{u1, u2} は正規直交系になる． 3. t3 = a3− su1− tu2 が{u1, u2} に直交するには 0 = (t3, u1) = (a3, u1)− s より s = (a3, u1), 0 = (t3, u2) = (a3, u2)− t より t = (a3, u2). ついで u3= t3/ p (t3, t3)とすれば{u1, u2, u3} は正 規直交系． 4. 以降，同様に続ける． 3DCGの場合 3DCGの場合は，式 (2) で求める。即ち，回転軸は求める直交行列の固有値 1 の固有ベク トルで，この固有ベクトルを含む一次独立な３個のベクトルを求めれば良い．2 番目のベクトルとして x 基 底ベクトルを第 1 候補として設定しておいて，回転軸が x-軸の場合は，y 基底ベクトルに置き換えればよ い．第 3 のベクトルは，第 1，第 2 のベクトルを直交化しておいて，それらの外積を取ればよい．以上で 回転軸を固有ベクトルとする正規直交系が求まるので，式 (2) を使用して回転行列が得られる．
1.3
空間の合同変換
空間で長さを変えない変換を合同変換 (congruent transformation) という．合同変換は直線を直線に写す ので線形変換となる．合同変換は平行移動 (translation)，回転変換 (rotation)，鏡映変換 (reflection) の合 成で得られる．1.4
その他の正則変換
次に，合同でない変換で，直線を直線に移す変換を考える．この場合も変換は行列で表されるが，長さは 保存されないので行列は直交行列ではない．行列の一般論から，正則性を仮定してそれを分類すると次の ようになる． 1.4.1 スケール変換 縦，横，高さの各方向の長さを変える変換をスケール変換 (scale transformation) という．縮尺比が同じ 場合は相似変換 (similar transformation) になる．行列の標準形を同次座標で表すと次の形をしている． S(a, b, c) = a 0 0 0 b 0 0 0 c (ただし abc6= 0) (3) 1.4.2 剪断変換1つの平面に平行な力をかけた変形をずれ変形または剪断変形 (shear or skew transformation) という．
剪断変形は不動平面に直交する座標軸を傾けることにより表現される．z–軸が剪断の軸の場合の行列を求 める．傾いた z–軸の単位ベクトルをt(α, β, 1)とすると，この写像は，次の行列で表される。 H = 1 0 α 0 1 β 0 0 1 −1 = 1 0 −α 0 1 −β 0 0 1 ここで，次の等式に注意する． 1 0 α 0 1 0 0 0 1 1 0 0 0 1 β 0 0 1 = 1 0 0 0 1 β 0 0 1 1 0 α 0 1 0 0 0 1 = 1 0 α 0 1 β 0 0 1 この式より剪断変換は x–軸，y–軸方向への剪断変換を任意の順に個別に行って合成すればよいことが分か る．特に，z–軸の傾斜角を用いてこの行列の成分を表すことにする． z–軸を x–軸正の方向へ θ，次に y–軸正の方向へ φ になるように剪断変換を行って，点t(0, 0, 1)がt(α, β, 1) に移ったとすると， α = tan θ β = tan φ であるので，この剪断変換の行列 H(θ, φ) は H(θ, φ) = 1 0 tan θ 0 1 tan φ 0 0 1 −1 = 1 0 − tan θ 0 1 − tan φ 0 0 1 (4)
となる．なお，z0–軸が x–軸，y–軸それぞれの負の方向へ傾いているときは，x–軸, y–軸の正の方向から
計った角度は 2π− θ であるが, これは tan(2π − θ) = tan(−θ) であるから, 座標軸の負の方向計った角を負
として計測すればよい. テキストによっては，z–軸の方からでなく，x–軸，y–軸の方から角度を測るように 指定したものがあるが，その場合は式 (4) の行列の要素の tan θ, tan φ をそれぞれ cot θ, cot φ に置き換え る．この場合, 角度の正負が分かりにくいので注意しよう. θ x x’ z tan θ z z’ x x’ z tan θ θ x z z’ z tan θ z tan θ x z z θ > 0 θ < 0 z z z’ y y’ z tan ϕ z tan ϕ y y’ y y z’ z ϕ ϕ z z tan ϕ z tan ϕ ϕ > 0 ϕ < 0
2
Hamilton
の四元数の利用
2.1
歴史と定義
複素数は実数の 2 次代数拡大体で，2 次元ベクトル拡大体でもある．3 次以上の多項式は可約であるか ら，3 次以上の代数拡大体はない．同じ理由で複素数の代数拡大体はない．次は 3 次元以上のベクトル拡 大体はないかと考えられる. 例えば，3 次元ベクトルのベクトル積は，同じベクトルのベクトル積は 0 ベク トルになるので，割り算が入らない．実際，次の定理により奇数次元の数は存在しない。 定理 単位元 e を持つ奇数次元の任意の (非可換) 体A は実数 R と同型であり、従って，1 次元である。証明 a∈ A を任意にとる。A における左積 La:A → A x 7→ ax はベクトル空間 A の線形写像である。 A は奇数空間であるから、写像 La は実数の固有値を持つ。その 1 つを λ、 対応する固有ベクトルを v6= 0 とおくと av = λv したがって、 (a− λe)v = 0. A は除法が可能で v 6= 0 あったから a = λe, すなわち、 a∈ Re. したがって A = Re. Hamiltonの 4 元数H は、3 次元の数はないかという問題の延長に発見された。Gauss も認識していたとい われている。しかし、残念なことに体の重要な性質である乗法の可換性は損なわれてしまう。そのために、 複素数のように日常的には使用されないのであるが、複素数の持つ幾何学的な性質を持っていて，これが computer graphicsにも利用可能である。特に，回転の計算において直交行列の計算より優位の面がある。 定義 4 次元実ベクトル空間 R4 の自然な基底 e0 = (1, 0, 0, 0), e1 = (0, 1, 0, 0), e2 = (0, 0, 1, 0), e3 = (0, 0, 0, 1) に次の乗法規則を定義する。 e0は単位元とする。 e1e1= e2e2= e3e3= e1e2e3=−e0
e1e2=−e2e1= e3, e2e3=−e3e2= e1, e3e1=−e1e3= e2
通常はベクトル記号を使用すると煩雑なので，e0= 1, e1= i, e2= j, e3= kを使用して α = a + ib + jc +
dk (a, b, c, d∈ R) と表す。
i2= j2= k2= ijk =−1 ij = −ji = k, jk = −kj = i, ki = −ik = j
とする。i, j, k に関する (複素数における虚数単位のような) 名称はないようである。 この記号を使用した積の公式は次の様になる。
(a + ib + jc + kd)(r + ix + jy + kz)
= (ar− bx − cy − dz) + i(ax + br + cz − dy) + j(ay − bz + cr + dx) + k(az + by − cx + dr)
2.2
他の表現方法
Hamiltonの四元数は上記のR4 の表現法の他に次の様に行列でも表現できる。 定理 次の写像はいずれも非可換体の準同型である． F : H → M(4, R) α = a + ib + jc + kd 7→ a −b −c −d b a −d c c d a −b d −c b a G : H → M(2, C) α = a + ib + jc + kd 7→ Ã a + ib −c − id c− id a− ib ! 注： 数の配列には幾つかのバリエーションがある。2.3
Hamilton
四元数の性質
Hamilton四元数の性質を列挙する。 (1) α = a + ib + jc + dk = a + q (ただし q = ib + jc + kd) の a を実部(またはスカラー部)、q を虚部(またはベクター部) と呼ぶ。a− q = a − ib − jc − kd を α の共役数と呼び、α と表す。 ImH = {ib + jc + kd|b, c, d ∈ R} を Hamilton 四元数の虚部空間と呼ぶ。 注： 虚部の定義が複素数の場合と異なる。 (2) α = a + ib + jc + kdに対し αα = a2+ b2+ c2+ d2. √αα =|α| と表し、|α| を α の絶対値と呼ぶ。 (3) α6= 0 のとき、上記より α−1= α |α|2. したがって、Hamilton 四元数は非可換体である。 (4) α (ただし|α| = 1) は次の一意的表現を持つ。 α = cos1 2ω + sin 1 2ωβ ただし β∈ ImH, |β| = 1, 0 < ω < 2π 証 α = a + q (ただし q ∈ ImH) とおく。b = |q|, β = q/b とおくと |β| = 1, α = a + bβ. αα = (a + bβ)(a− bβ) = a2+ (bβ)2= a2+ b2= 1より a = cos1 2ω, b = sin 1 2ω をみたす 0 < ω < 2π が唯一つ存在する。 注： この式は複素数の極形式 (Euler from) に相当する。以後、この形式を四元数の標準形式と呼ぶ。 一方，i, j, k を使用した形式は自然な形式と呼ぶことにする。2.4
幾何学への応用
四元数を用いて、3 次元ユークリッド空間の回転を扱う。始めに次の式を考える。 α = a + ib + jc + kd (ただし|α| = 1) ζ = w + ix + jy + kzαζα = w + i{(a2+ b2− c2− d2)x + 2(bc− ad)y + 2(ac + bd)z} + j{2(bc + ad)x + (a2− b2+ c2− d2)y + 2(cd− ab)z} + k{2(bd − ac)x + 2(ab + cd)y + (a2− b2− c2+ d2)z}
αζ α = w + ix0+ jy0+ kz0 と置くと，この式は次の様に行列で表示できる。 x0 y0 z0 w = a2+ b2− c2− d2 2(bc− ad) 2(ac + bd) 0 2(bc + ad) a2− b2+ c2− d2 2(cd− ab) 0 2(bd− ac) 2(ab + cd) a2− b2− c2+ d2 0 0 0 0 1 x y z w (1) この行列は列ベクトルが正規直交系を作るので，直交行列である。
αを標準形式 α = cos1 2ω + sin 1 2ω β, β = ip + jq + kr, p 2+ q2+ r2= 1 にして、上記の行列の左上 3× 3 部分を T とおくと次の様になる。ただし ω0= ω 2 とおく。 T =
cos2ω0+ (2p2− 1) sin2ω0 −2r cos ω0sin ω0+ 2pq sin2ω0 2q cos ω0sin ω0+ 2pr sin2ω0
2r cos ω0sin ω0+ 2pq sin2ω0 cos2ω0+ (2q2− 1) sin2ω0 −2p cos ω0sin ω0+ 2qr sin2ω0 −2q cos ω0sin ω0+ 2pr sin2ω0 2p cos ω0sin ω0+ 2qr sin2ω0 cos2ω0+ (2r2− 1) sin2ω0
T について次のことがいえる。 (2) tT T = E (3) t(p, q, r)が固有値 1 の固有ベクトルである。 (4) Tr(T ) = 3 cos2ω0− sin2ω0 = 1 + 2 cos ω
である。 これらのことから、行列 T がt(p, q, r)が回転軸、回転角が ω の回転行列であることを表している。最後 のことは T の固有方程式の 2 次の係数が Tr(T ) = 1 + eiω+ e−iω であることによる。
2.5
CG
への応用
これまで説明したことから回転の定義に四元数が使えることが分かったが，通常の方法に比べ優位であ る点があるかを考えてみる。始めに，四元数の計算を直接行うプログラムは存在しないから、四元数の表 示形式が簡単であることは、計算が高速になる理由にはならない。通常の方法では、回転軸 ~p1 と回転角 θ を指定されたときに、次の行列を計算する必要がある。 直交行列 P = (~p1, ~p2, ~p3), S = 1 0 0 0 cos θ − sin θ 0 sin θ cos θ とし、P SP−1. ただし、~p1, ~p2, ~p3は正規直交系。このときに、~p2, ~p3は ~p1 から直接計算する方法がなく，(1) 適当なベク トル，例えば，座標軸ベクトル，を 2 本追加してシュミットの直交化法で求める，(2) 適当なベクトル，例 えば，座標軸ベクトル，とベクトル積を取り，正規化して ~p2，さらに ~p1, ~p2 のベクトル積で ~p3を求める． (3) ~p1 に直交するベクトルを連立方程式で求める，などの方法が考えられるが，いずれも曖昧性が発生し， その排除の判定が必要になるので，計算が遅くなる． 一方，四元数の場合は α = a + ib + jc + kd を指定した場合に、β = ib + jc + kd とおいて、α = a + β, cos θ/2 = a/|α|, sin θ/2 = |β|/|α| により回転角 θ、t(b, c, d)により回転軸が指定されたことになり、この 定数を用いて上記の行列 T が曖昧性なく決定される。従って，正規直交系を決定する時間だけ高速になる ことになり、これがアニメーションのように行列を連続的に求めなければならないときに smoothness の点 で優位になると評価されているのであろう。 もっとも、四元数により回転を指定する場合は、予め使用者の方で計算しておかなければならない部分が 増える。2.6
メモ
以上の説明は、参考書を見て再構成したものであるが、それに至るまでの経過も何かの参考になると思わ れるので、メモにしておく。
始めに、四元数が CG で使用されていることは、例えば，“Interactive Computer Graphics, E. Angel” に上述の結果が記載されている。一般には、計算は行列法より高速になると説明されている．アフィン 空間と比較して虚部空間がユークリッド空間に対応していると見て、四元数 2 つの積を見てみる。単純 積だけでは実部が変わってしまうので、αβα で打ち消して考えればよいとする。α = i cos θ + j sin θ と
して計算してみると、回転軸がt(cos θ, sin θ, 0) であることが得られたが，回転角が π になった。回転
角は指定してないと思ったので，実部が影響していると考えて、逆に回転角を指定して対応する四元数 をもとめると、α = cos ω/2 + i sin ω/2 cos θ + j sin ω/2 sin θ が得られた。ここまでの結果を見て、α = cos ω/2 + (ip + jq + kr) sin ω/2、(ただし、p2+ q2+ r2= 1、即ち，これが標準形式であることを知らな
かったのであるが)、を用いて上記の結果が得られた。多項式環の参考書を調べて、改めて標準形式の用語 と結果の正しさを確認して、上記の様に纏めることができた。
3
ユークリッド空間とアフィン空間
3.1
はじめに
CGでは平行移動，回転，剪断，投影，射影など直線を直線に移す様々な変換が扱われる。これらの変換 は，行列の積，ベクトルの和，一次分数変換等で表現されるので，別々の関数とパラメータを用意しておか なければならない。しかし，ユークリッド空間をアフィン空間に埋め込めば，これらの変換は同次座標を用 いたアフィン写像という行列の掛け算だけで達成できることをことを説明し， CG の構成に必要な計算を 求める。3.2
アフィン空間
ベクトル空間 V と集合 A において，任意の a∈ V と p ∈ A の間に和 p + a ∈ A が定義されて，次の条件が成立するとき，A を V を基準ベクトル空間 (standard vector space) とするアフィン空間 (affine
space)であるという． (1) (p + a) + b = p + (a + b) (p∈ A, a, b ∈ V) (2) 任意の p, q∈ A に対し，q = p + a となる a ∈ V が唯一つ存在する． 例 (i) ベクトル空間はそれ自身を基準ベクトル空間とするアフィン空間である． (ii) n + 1次元のベクトル空間 W に対し，n 次元の部分空間 V を取り，A ={p + a | 0 6= p 6∈ V, a ∈ V} とおくと A は V を基準ベクトル空間とするアフィン空間である． この例では，o∈ A を任意に固定すると p ∈ A は p = o + a と一意に表されるので，a を o を始点
V を基準ベクトル空間とするアフィン空間 A の r + 1 個の点 pi (0 5 i 5 r) に対し，r 個のベクトル −−→ p0pi (15 i 5 r) が V で線形独立であるとき，これらの点は独立であるといい，そうでないとき従属であ るという．基準ベクトル空間の次元をそのアフィン空間の次元という．アフィン空間の次元が n であると き，独立な点の最大個数は n + 1 である． Vを基準ベクトル空間とする n 次元アフィン空間があるとき，e1, . . . , en を V の一つの基底，o∈ A を 任意に固定元とすると，任意の p∈ A は p = o + n X i=1 xiei と一意に表される．このF = (o ; e1, . . . , en)を A のアフィン枠 (affine frame) という． これをベクトル形式で表すときは次のように表し，同次座標 (homogeneous coordinate) ともいう． x1 .. . xn 1 A, B を V, W をそれぞれ基準ベクトル空間とするアフィン空間とする．写像 ϕ : A→ B がアフィン写像 (affine mapping)であるとは，ϕ に対し，線形写像 ψ : V→ W があって ϕ(p + x) = ϕ(p) + ψ(x) がすべ ての p∈ A, x ∈ V に対し成立するときをいう．さらに，A = B で ϕ が全単射のときは ϕ をアフィン変換 (affine transformation)という．アフィン写像は同次座標系を用いて次の形式で表される． Ã M t 0 1 ! ここで前説で説明したユークリッド空間の合同変換に対応するアフィン変換行列を整理しておく。 (1) 直交変換，スケール変換，剪断変換 P に対し Ã P t0 0 1 ! (2) 平行移動 t = (α, β, γ) に対し Ã E tt 0 1 ! ただし 0 = (0, 0, 0)，E は単位行列である。
3.3
アフィン枠の変換に伴う座標変換
Aを V を基準ベクトル空間とするアフィン空間，(o ; e1, . . . , en), (o0 ; f1, . . . , fn)をアフィン枠とする． アフィン枠の変更に伴う座標の変換を求める．p∈ A を任意にとり，それぞれのアフィン枠を使うと p = o + n X i=1 xiei = o0+ n X j=1 yjfj o0 = o + n X i=1 aiei ここで， fj = n X i=1 tjiei (j = 1, . . . , n) とおいて以上の式を代入すると p = à o + n X i=1 aiei ! + n X j=1 yj à n X i=1 tjiei ! = o + n X i=1 ai+ n X j=1 yjtji ei 従って xi= ai+ n X j=1 yjtji (i = 1, . . . , n) 同次座標を用いて行列の形式で表すと x1 .. . xn 1 = t11 . . . tn1 a1 .. . . .. ... ... t1n . . . tnn an 0 . . . 0 1 y1 .. . yn 1 または y1 .. . yn 1 = t11 . . . tn1 a1 .. . . .. ... ... t1n . . . tnn an 0 . . . 0 1 −1 x1 .. . xn 1 (5)
3.4
3DCG
の場合
3DCGの場合は，物体や変換を記述する座標系をワールド枠という。ワールド枠で注意しなければいけ ないことは，CG の場合の各座標軸の取り方で，x–軸は横軸，y–軸は縦軸，z–軸は前後軸になり，右，上，前が夫々正の方向である。数学では縦軸は z–軸である。一方，投影の計算をワールド枠で行うと投影の理 論の説明が困難になる。そのために，3DCG ではカメラ枠という投影用のアフィン枠を別に用意し，モデ リングはワールド枠で，レンダリングはカメラ枠で行い，両者の変換をアフィン変換で行う。
カメラ枠は，原点 VRP(View–Reference Point), 投影面の縦軸 VUP(View–Up Vector), 投影面の法線
VPN(View–Plane Normal)で定義される。VUP がカメラ枠の y–軸，VPN が z–軸になる。ただし，VPN
の方向は負の方向とする。x–軸はこの 2 つのベクトルのベクトル積を取ればよい。この座標軸の取り方は ワールド枠の座標軸の取り方と同様で，縦軸が y–軸になっている。これらをワールド枠で定義したときに， VRPが (ex, ey, ez)，VUP が (ux, uy.uz)，VPN が (nx, ny, nz)，(α, β, γ) = (ux, uy, uz)× (−nx,−ny− ny) であるとき，アフィン変換行列 T を T = α ux −nx ex β uy −ny ey γ uz −nz ez 0 0 0 1 とおくと，ワールド枠で座標が (wx, wy, wz)の点をカメラ枠の座標で (cx, cy, cz)であるとすると，上記の 説明から T cx cy cz 1 = wx wy wz 1 または cx cy cz 1 = T −1 wx wy wz 1 の関係が成り立つ。 OpenGLの場合，この変換を行う関数は
gluLookAt(eyex eyey, eyez, atx, aty, atz, upz, upy, upz)
であって，(eyex, eyey, eyez) が視点 (COP)，(atx, aty, atz) が着目点，即ち，レンズの軸方向にある 空間の点，(upx, upy, upz) が上方ベクトルであるから VRP = (eyex, eyey, eyez), VUP = (upx, upy, upz), VPN = (atx− eyex, aty − eyey, atz − eyez) で VUP と VPN をこの順に取ったベクトル積を (α, β, γ) と すると T =
α upx eyex− atx eyex β upy eyey− aty eyey γ upz eyez− atz eyez
0 0 0 1 である
3.5
OpenGL
における行列についての注意
3.5.1 行列積の順序 OpenGLでは, 図形の変形の場合でも投影の変換の場合でも同次座標とアフィン写像で行われる．OpenGL では. 回転，平行移動，スケール変換，剪断変換を定義する場合に，行列を掛ける場合に postmultiplication で行われる．すなわち，プログラム文中で行列を記述すると，現在の行列 (CTM = current transformation matrix)の後ろに記述された行列を掛けたものが新しい CTM になる．実際の変換は x0 = Mx (M : CTM) であるから，2 つ以上の行列を掛けて変換を定義する場合は記述する順序に注意しなければならない．例え ば，原点以外の固定点 a の周りの回転を実施する場合は，固定点から原点への平行移動の行列 T(−a)，原 点の周りの回転の行列 R(θ)，原点から固定点への平行移動の行列 T(x)，の順に変換を行う．式で表すと x0 = T(a)R(θ)T(−a)x であるが，プログラム文では glMatrixMode(GL_MODELVIEW); glLoadIdentity();glTranslatef(ax, ay, az); glRotatef(t, rx, ry, rz); glTranslatef(-ax, -ay, -az); /* define objects here */ glutSwapBuffer(); となり，式とプログラム文の順序は同じで，言葉での説明と逆であることに注意． 3.5.2 行列の定義と導入法 OpenGLで扱う行列はすべて 4× 4 行列であるが，API が持っていない行列をユーザが新たに定義する には 2 次元配列により行列を定義するのでなく，1 次元配列で定義しなければならない．1 次元配列から行 列への変換は，先頭から順に第 1 縦 (列) ベクトル，第 2 ベクトル，というように代入される．配列の指数 を行列の形で表すと 0 4 8 12 1 5 9 13 2 6 10 14 3 7 11 15 例えば， M = 1 2 0 0 0 1 3 0 0 0 1 0 0 0 0 1
の場合は， GLfloat m[16]; for(i=0;i<15;i++) m[i]=0.0; m[0]=m[5]=m[10]=m[15]=1.0; m[4]=2.0; m[9]=3.0; となる． ユーザが定義した行列を使って計算したい場合は，次のいずれかのプログラム文を使用する．新しい行列 を CTM に定義したい場合は glLoadMatrixf(myarray) CTMにユーザ定義の行列を掛けたい場合は glMultMatrixf(myarray) また，計算途中で CTM を一時変更して計算終了後元に戻したい場合は次のようにスタックを利用する． glPushMatrix(); glTranslatef(...); glRotatef(...); glScalef(...); /* draw objects */ glPopMatrix();
4
CG
と投影
(Projection)
4.1
投影モデル
CG における造影はカメラによる造影をモデルとしている．レンズを投影の中心 (COP = center ofprojection)として，物体と COP を結ぶ直線を投影線 (projector) といい，造影は投影線とフィルムの交点
を物体の投影点として行われる．即ち，ピンホールカメラによる造影をモデルにして考える．このままでは 作られる像が上下，左右が反転するので，フィルム面はレンズの前にあるものとして投影点を求める．さら に，実際のカメラと異なり，コンピュータによる処理は有限の範囲にしなければならないので，投影の左 右，上下，前後枠を設定する．この枠を視野領域 (view volume) またはクリッピング (clipping) という．カ メラで撮影する感覚では左右と上下の制限は自然であるが，前後の制限特に前面の感覚はないので注意が 必要である．前後のクリッピングは芸術的な写真では「ボケ」などといい，重要な要素ではあるが，CG で は視野領域内にあるものはボケがなく写されて，逆にボケの状態を作る方が難しい。
4.2
投影の種類
4.2.1 平行投影 (parallel projection)
投影線が互いに平行になっている場合である。これは COP が無限遠にある場合という説明もできるが，
COPの代わりに投影直線の方向ベクトル DOP(Direction Of Projection) をパラメータに使う。この平行
投影には次の 2 通りがある。
直投影，正射影 (orthogonal or orthographic projection) 投影線と投影面が直交する場合で，DOP
が投影面の法線ベクトルと平行になる。この投影では投影面と物体との距離に関係なく同じサイズに投影 される．
斜め投影 (oblique projection) 投影線と投影面が直交しない場合で，DOP は投影面の法線ベクトルと
は平行にならない。この投影では平行性は保たれるが，サイズは方向により保存されない。蛇腹式 (bellows) カメラでフィルム面をレンズに対し斜めにした場合に相当する． 4.2.2 透視投影 (perspective projection) 物体と COP の中間に投影面があるものとし，投影線とフィルム面の交点を投影点とする．この投影では， 直線は直線に投影されるが，平行性は保存されない．同じサイズの物体でも投影面に近いものが遠くにあ るものより大きく投影される．即ち, 遠近法の表現が得られることになり, これにより立体感が表される． 人の手で 3 次元の図を描くときの透視投影法では, 平行線の表現方法に 3 種類の描画法が考えられる．1 番目は 1 点透視法で，例えば，建物の真正面に立って見た感じを表す方法で，縦横の平行線は平行に描き， 前後方向の平行線は正面の中心に収束するように描く方法である．2 番目は 2 点透視法で，建物の角に立っ て見た感じを表す方法で，縦の平行線は平行に描かれるが，横の平行線は中心より左側の平行線は左の 1 点 に，右側の平行線は右の 1 点に収束するように描く方法である．3 番目は通常の透視法で，建物の例では角 の空中から見た感じを表す方法で，横方向の平行線は左側の平行線は左側の 1 点に，右側の平行線は右側 の 1 点に，縦方向の平行線は，下から見上げた場合は上，上から見下ろした (俯瞰) 場合は下の１点に収束 するように描く方法である．CG でこれらの違いを表すには，カメラ枠の設定で表現可能である．
5
投影の計算と正規化
5.1
作像の手順
OpenGLでは，作像を物体の変換 (tranform) と投影 (projection) の 2 つのプロセスに分けて行い，それぞ
れの変換に必要な行列は，それぞれの行列スタックに定義される。具体的には，変換からカメラ枠までの行列 は glMatricMode(GL_MODELVIEW) にスタックされ，投影に関する行列は glMatrixMode(GL_PROJECTION) にスタックされる。これらのプロセス間ではデータは行ったり来たりすることなく，最初のプロセスで得ら れたデータは次のプロセスに連続して渡されて処理されて投影像を作成する。このようなプロセスの連結 方法をパイプラインと呼んでいる。これらのプロセスは実際はアフィン枠を使用したことにより同次座標 に対する変換行列の積とユークリッド座標への引き戻し (代表系の取り出し) で実現される．
5.2
投影
投影の計算を求める．計算はすべてカメラ枠を用いて行われる．
5.2.1 平行投影
平行投影は投影線が互いに平行である投影である．COP に代って DOP がパラメータに使われる。
1. 直投影 直投影，正射影 orthogonal projection または orthographic projection は投影線が投影面 (フィ ルム面) に直交している．DOP は投影面に直交している，即ち，カメラ枠の z–軸に一致している。直投影 の場合は投影はフィルムの位置に依存しないので，フィルム面の方程式は z = 0 としてよい．したがって， 投影変換は次の行列になる． xp yp zp 1 = Porth x y z 1 = 1 0 0 0 0 1 0 0 0 0 0 0 0 0 0 1 x y z 1 実際は (zp= 0) 2. 斜め投影 斜め投影は DOP が投影面に直交していないが，決して特殊な投影ではない．例えば，3 次 元の図を描くときに数学の通常の表現では，横軸を y–軸，縦軸を z–軸とし直交させ，x–軸を斜め手前に取 る．CG では x, y–軸が横縦の軸で，z–軸を斜めに取るが，これは斜め投影をしていることになる． OpenGLでは斜め投影を直接定義する関数がないので，投影面に平行な方向へ斜め投影をしてから直投影 をするプロセスを取る．斜め投影の部分はユーザが計算する必要がある．斜め投影は z–軸を剪断軸とする 剪断変換と考えることができる．斜め投影を DOP で特徴づける。DOP の xz–平面への投影と z–軸の成す 角を θ，yz–平面への投影と z–軸の成す角を φ，とする．角度の正負については, 剪断変換の場合と同じよう に考える. 即ち, もし，DOP が x–軸，y–軸の正の方向へ傾いている場合は，角度は正の角であるが，負の 方向へ傾いている場合は，負の角と定める．この投影で (x, y, z) の投影点を (xp, yp, zp) (ただし zp= 0) とすると xp− x zp− z = tan θ したがって xp= x− z tan θ 同じように yp− y zp− z = tan φ したがって yp= y− z tan φ これより次のように斜め投影の行列を求めることができる．これは式 (4) と同じである． H = H(θ, φ) = 1 0 − tan θ 0 0 1 − tan φ 0 0 0 1 0 0 0 0 1
剪断変換のときに説明したように，テキストによっては角度を DOP と z–軸の成す角でなく，DOP と x– 軸，y–軸の成す角度を用いたものがある．その場合は，行列の成分中の tan θ, tan φ をそれぞれ cot θ, cot φ と換えればよい．その場合，角度の符号に関しては同じ注意を守ることに注意しよう． 以上より，斜め投影の行列は次のようになる． P = PorthH 具体的には，matrixH を H に対応するものとして，次のようになる. glMatrixMode(GL_PROJECTION); glLoadIdentity(); glOrtho(-4.0,6.0,-4.0,6.0,-20.0,20.0); glMatrixMode(GL_MODELVIEW); glLoadIdentiry(); gluLookAt(0.0,0.0,6.0,0.0,0.0,0.0,0.0,1.0,0.0); glMultMatrixf(matrixH);
/* define view objects */
glOrthoのパラメーターが対象でない点にも注意しよう. 斜め投影の場合の視点は z–軸上にないと斜め投 影の特徴が消えてしまう。 x z DOP θ oblique skew ortho y z DOP oblique ortho skew ϕ θ > 0 ϕ < 0 5.2.2 透視投影 透視投影の原理は COP と投影面を設定し，物体と COP を結ぶ直線 (投影線) と投影面の交点を求める ことである．規定の位置は，COP が原点，投影面は z = d (d < 0) である．透視投影では平行投影と異な り，view volume は COP と頂点とする四角錐を切り取った台 (frustum) になる．
カメラ枠を使用して，点 (x, y, z) を (xp, yp, zp)に透視投影する場合の計算式は次のようになる．まず， 投影面の方程式より zp= d 次に (x, y, z) と (xp, yp, zp)を結ぶ直線が原点を通ることから，各座標は比例しているので， xp x = yp y = zp z
したがって， xp= x z/d yp= y z/d z = d ここで，ベクトル表示すると xp yp zp = x z/d y z/d d であるから，分数変換を含み行列の計算ではこの結果を得ることができない．しかし，同次座標を使用す ると 1 0 0 0 0 1 0 0 0 0 1 0 0 0 1/d 0 x y z 1 = x y z z/d → 第 4 成分による割り算 → x z/d y z/d d 1 と計算でき，最後の代表元を取り出す計算は z 成分には独立の計算であるから予め組み込むことが可能で ある．これで同次座標を利用するする利点が理解できるであろう． OpenGL における透視投影は次のいずれかの関数を使用する． glFrustum(left,right,bottom,top,near,far) gluPerspective(fovy,aspect,near,far)
gluPerspectiveにおける f ovy は field of view または angle of view の意味で，clipping の縦軸を底辺と する三角形の頂点 COP の頂角であり，aspect は clippping の横対縦の比 w/h である．near, f ar の意味 は glFrustum に同じである．
5.3
Normalization
5.3.1 平行投影の場合 OpenGL では，直投影の関数は次のように宣言される． glMatrixMode(GL_PROJECTION); glLoadIdentity(); glOrtho(left,right,bottom,top,near,far);この宣言で near と far はカメラレンズからの距離を意味するので，0 < near < f ar でなければならない． この宣言では clipping が次のようになる．
lef t≤ x ≤ right, bottom ≤ y ≤ top, −near ≤ z ≤ −far
defaultの clipping は
なので，上で宣言された clipping が表示面に正しく表示されるために，(フィルム面に目的の世界が収まる
ようにするために)，フィルム枠が x =±1, y = ±1, z = ±1 となるように座標を変換する必要がある．こ
れを正規化 (normalization) という．正規化の行列は，まず clipping の中心をカメラ枠の中心にするため に (−(xmax+ xmin)/2,−(ymax+ ymin)/2,−(zmax+ zmin)/2)平行移動する．次に，clipping が
x =±1, y =±1, z =±1
になるように x–軸方向で 2/(xmax− xmin)倍，y–軸方向で 2/(ymax− ymin)倍，z–軸方向で 2/(zmax− zmin)
倍する．したがって，求める変換は P = PorthST = Porth 2 xmax−xmin 0 0 0 0 2
ymax−ymin 0 0
0 0 2 zmax−zmin 0 0 0 0 1 1 0 0 −xmax−xmin 2
0 1 0 −ymax−ymin
2 0 0 1 −zmax−zmin 2 0 0 0 1 = Porth 2 xmax−xmin 0 0 − xmax+xmin xmax−xmin 0 y 2
max−ymin 0 −
ymax+ymin
ymax−ymin
0 0 2 zmax−zmin − zmax+zmin zmax−zmin 0 0 0 1 上のパラメータはカメラから見たもので，カメラレンズが z–軸の負の方向が正面であることを考慮すると， このパラメータを使用するときは左手系に変換しなければならない．即ち，z 成分を−1 倍する必要があ
る．上記のアフィン行列は zmax=−near, zmin=−far を代入して z 成分を −1 倍して次を得る．
P = Porth 2 right−left 0 0 − right+lef t right−left 0 2 top−bottom 0 − top+bottom top−bottom 0 0 − 2 f ar−near f ar+near f ar−near 0 0 0 1 5.3.2 Perspective関数の場合 gluPerspective(fovy,aspect,near,far) 関数を使用して透視投影を行う際の正規化について考える．
透視投影の規定値は，f ovy = 90◦ で view volume の側面が投影線と 45◦をなすように設定されている．言
い換えれば，view volume は x =±z, y = ±z で z 負の方向が無限のピラミッド型の領域である．そのた めに，near, f ar のパラメータで z 方向の制限を設けてあるが，これも z =±1 に正規化したい． まず，投影面を z =−1 とすると，上で求めた投影の行列で d = −1 とおくと M = 1 0 0 0 0 1 0 0 0 0 1 0 0 0 −1 0
ここで次の行列を考える． N = 1 0 0 0 0 1 0 0 0 0 α β 0 0 −1 0 αβ6= 0 とするとこの行列は正則である．この行列 N と直投影の行列 Morthを合成すると MorthN = 1 0 0 0 0 1 0 0 0 0 0 0 0 0 −1 0 この変換を用いると x0 y0 z0 w0 = MorthN x y z w = x y 0 −z 最後に第 4 成分による割り算 (perspective division) を経て xp=− x z yp=− y z という透視投影を得る．一方， x0 y0 z0 w0 = N x y z w より x0= x, y0= y, z0= αz + β, w0 =−z であるから perspective division を行うと x00=−x z, y 00=−y z, z 00=− µ α +β z ¶ ここで − µ α + β zmax ¶ = 1 − µ α + β zmin ¶ =−1 すなわち α =−zmax+ zmin zmax− zmin β = 2zmaxzmin zmax− zmin
とおけば，z =±1 という正規化が得られる．この値を代入した行列 N を透視投影の正規化行列という．こ のとき z00=− µ α +β z ¶ という変換は分数変換であるが， zmax, zmin ともに z = 0 に対して同じ側に設定されるから β > 0 で，し たがって z1< z2 =⇒ z001 < z002 の関係が成立するので，割り算という計算精度上の問題があるが，この変換で z–バッファの機能は正常に 働くことが保証される．唯一の問題点は，この変換により始めの投影面 z =−1 が z = β − α に変わるこ とであるが，実際の投影はこの後に直投影 Morth を施すことから，投影図の変形はおきない． 5.3.3 Frustum関数の場合 glFrustumの規定値は x =±z, y = ±z, −1 ≤ z ≤ 1 であり，glFrustum(left,right,bottom,top,near,far) のパラメータをそのまま使用すると，ディスプ レイの隅の方にのみ画像が現れることになる．悪くすると画像が現れない．これに対し，カメラをそちらの 方へ向け直すとフィルム面が移動して造影法が異なることになる．この場合の正規化は，入力されたパラ メータの中心を通る軸が z–軸になるように剪断変換を行い，次にスケール変換を行う．最後に，z 方向の 正規化を行う． まず，上記のパラメータによる投影の中心は µ xmin+ xmax 2 , ymin+ ymax 2 , zmax ¶ であるから傾斜した軸の方向ベクトルはt((x
max+ xmin)/2zmax, (ymax+ ymin)/2zmax, 1)であり，したがっ
て最初の剪断変換は，式 (4) より H = 1 0 xmax+xmin 2zmax 0 0 1 ymax+ymin 2zmax 0 0 0 1 0 0 0 0 1 −1 = 1 0 −xmax+xmin 2zmax 0 0 1 −ymax+ymin 2zmax 0 0 0 1 0 0 0 0 1 θ, φを次のようにする．
cot θ =−xmin+ xmax 2zmax
, cot φ =−ymin+ ymax 2zmax 上記の剪断変換の行列は H = H(θ, φ) = 1 0 cot θ 0 0 1 cot φ 0 0 0 1 0 0 0 0 1
この座標変換で透視投影の view volume の COP に近い窓が H(θ, φ) xmin ymin zmax 1 = xmin−xmax 2
ymin−ymax
2 zmax 1 H(θ, φ) xmax ymax zmax 1 = xmax−xmin 2
ymax−ymin
2 zmax 1 即ち，frustum が x =±xmax− xmin 2zmax z y =±ymax− ymin 2zmax z zmin≤ z ≤ zmax である．これを x =±z y =±z とするために，次のスケール変換を行う． S( 2zmax xmax− xmin , 2zmax xmax− xmin , 1) 最後に z =±1 にするための正規化 N を行う．これを以上の順に合成すると，次の投影行列 P を得る．なお，この変換 は API が自動的に行うので，ユーザが扱う必要はない． P = NSH = 2zmax xmax−xmin 0 − xmax+xmin xmax−xmin 0 0 2zmax
ymax−ymin −
ymax+ymin
ymax−ymin 0
0 0 −zmax+zmin zmax−zmin 2zmaxzmin zmax−zmin 0 0 −1 0 さらに OpenGL のパラメタを用いてこの正規化行列を求めると
zmax=−near, zmin=−far
を代入して， OpenGL では z–軸方向がカメラ枠の逆向きであったから，第 3 行を−1 倍して次の行列を 得る． P = 2zmax xmax−xmin 0 − xmax+xmin xmax−xmin 0 0 2zmax
ymax−ymin −
ymax+ymin
ymax−ymin 0 0 0 −f ar+nearf ar−near −2f arf ar−near×near
0 0 −1 0
6
平面への影を求める
6.1
平行光線による影
平行光線の方向ベクトルを (p, q, r) ，影のできる平面の方程式を lx + my + nz + c = 0 とする．点 P = (x, y, z) の上の面の影の座標を P0= (xs, ys, zs)とすると，P0 が影の面上にあることより lxs+ mys+ nzs+ c = 0 −−→ P0Pが光源の方向ベクトルに一致することより xs− x p = ys− y q = zs− z r 最後の式を k とおいて，影平面の方程式へ代入して k =−lx + my + nz + c K0 ただし K 0 = lp + mq + nr したがって xs = (K0− lp)x − mpy − npz − cp K0 ys = −lqx + (K 0− mq)y − nqz − cq K0 zs = −lrx − mry + (K0− nr)z − cr K0 これを同次座標を用いて次のように行列の計算で表される． Moshd x y z 1 = K0− lp −mp −np −cp −lq K0− mq −nq −cq −lr −mr K0− nr −cr 0 0 0 K0 x y z 1 (6) これを次のように分解する． Moshd= K0E4− p 0 0 0 0 q 0 0 0 0 r 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0 l 0 0 0 0 m 0 0 0 0 n 0 0 0 0 c (7) しかし，OpenGL には行列の加減を行う関数がないので，式 (6) を直接代入するか，式 (7) の計算の関数 を別に用意する必要がある．特に影のできる平面が y + c = 0 で上からの平行光線の場合は l = n = p = r = 0, m = 1, q =−1であるから，同次座標の代表を取り出す計算を考慮して Moshd x y z 1 = 1 0 0 0 0 0 0 −c 0 0 1 0 0 0 0 1 x y z 1 = x −c z 1 また，影平面が y + c = 0 で (p, q, r) 方向からの平行光線の場合は l = n = 0, m = 1 であるから，同様に M x y z 1 = q −p 0 −cp 0 0 0 −cq 0 −r q −cr 0 0 0 q x y z 1 となる.
6.2
点光源による影
点光源の位置を L = (p, q, r)，影のできる平面の方程式を lx + my + nz + c = 0 とする．点 P = (x, y, z) の影平面上の点を P0= (xs, ys, zs)とすると P0 が影平面上にあることから lxs+ myx+ nzs+ c = 0 L, P, P0 が 1 直線上にあることから (xs− p, ys− q, zs− r) = k(x − p, y − q, z − r) これを影平面の式に代入して k =− K lx + my + nz + c− K ただし K = lp + mq + nz + c これより xs = (K− lp)x − mpy − npz − cp K− (lx + my + nz + c) ys = −lqx + (K − mq)y − nqz − cq K− (lx + my + nz + c) zs = −lrx − mry + (K − nr)z − cr K = (lx + my + nz + c)したがって影を作る行列は Mpshd x y z 1 = K− lp −mp −np −cp −lq K− mq −nq −cq −lr −mr K− nr −cr −l −m −n K− c x y z 1 これを次のように分解する． Mpshd= KE4− p 0 0 0 0 q 0 0 0 0 r 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 l 0 0 0 0 m 0 0 0 0 n 0 0 0 0 c (8) この行列もこのまま利用すると，OpenGL では行列の計算を行う関数を別に用意しなければならないが，次 の変形を考える． 1 0 0 −p 0 1 0 −q 0 0 1 −r 0 0 0 1 p 0 0 0 0 q 0 0 0 0 r 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 l 0 0 0 0 m 0 0 0 0 n 0 0 0 0 c 1 0 0 p 0 1 0 q 0 0 1 r 0 0 0 1 = p 0 0 −p 0 q 0 −q 0 0 r −r 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 l 0 0 lp 0 m 0 mq 0 0 n nr 0 0 0 c = 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 l 0 0 lp 0 m 0 mq 0 0 n nr 0 0 0 c = 0 0 0 0 0 0 0 0 0 0 0 0 l m n K したがって 1 0 0 −p 0 1 0 −q 0 0 1 −r 0 0 0 1 Mpshd 1 0 0 p 0 1 0 q 0 0 1 r 0 0 0 1 = K 0 0 0 0 K 0 0 0 0 K 0 −l −m −n 0 したがって Mpshd= 1 0 0 p 0 1 0 q 0 0 1 r 0 0 0 1 K 0 0 0 0 K 0 0 0 0 K 0 −l −m −n 0 1 0 0 −p 0 1 0 −q 0 0 1 −r 0 0 0 1 (9) この形にすれば OpenGL の行列の掛け算で点光源による影が簡単に求めることができる．実は，式 (9) の 中央の行列は透視投影の行列であり，点光源による影は影を投影する立地を光源毎，光源が原点にくる様に
平行移動して，原点からの透視投影で影を求めて，元の位置に戻して求めれば計算が簡単になることを示 している。平行光線による造影も同じ考えで得られるが，この場合は真上からの光線以外のときは斜め投 影に当たるので，行列は上記の様に単純ではない。 特に，影の平面を y + c = 0 即ち l = n = 0, m = 1 とすれば，行列 (9) は Mpshd= 1 0 0 p 0 1 0 q 0 0 1 r 0 0 0 1 1 0 0 0 0 1 0 0 0 0 1 0 0 −q+c1 0 0 1 0 0 −p 0 1 0 −q 0 0 1 −r 0 0 0 1 (10)
6.3
造影のプログラム
影を作るプログラムは，上記で説明した様に定義することができる。一般には，点光源による造影が光源 に関する平行移動を組み合わせた方が分かり易いという説明がなされることが多いが，平行光線による造 影行列 (7) と点光源による造影行列 (8) は，第 4 行だけの違いであり，行列 (9) を用いるにしても行列の計 算回数は，行列 (8) の場合と同じである。従って，造影プロセスを導入する場合に，行列の掛け算プログラ ムを 1 つ用意しておけば，平行光線，点光源いずれの場合も同じ考えで造影ができる。 次の例では，行列 (10) により点光源による影を作成するもので，本体は transform も含め，独立のプロ シージャ object() で定義しておく。光源が固定されている場合は，影を投影する行列 shdm[] は main() の先頭部分など描画ループを考慮した場所におく。 void o_shadow() { GLfloat shdm[16]; for(i=0;i<15;i++) shdm[i]=0; shdm[0]=shdmm[5]=shdm[10]=1; shdm[7]=-1/(q+c); glDisable(GL_LIGHT0); glPushMatrix(); glTranslatef(p,q,r); glMultMatrix(shdm); glTranslatef(-p,-q,-r); object(); glPopMatrix(); glEnable(GL_LIGHT0); } この方法で得られる影は，影の境界がくっきりしていて，色も真っ黒であり，現実の影とは違って見える。 色については，glDisable(GL_LIGHT0) を止めて，影も通常のオブジェクトと同様に材質光を指定すれば よいが，このときに本体オブジェクトの材質光と同等の扱いをする必要がある。具体的には，影を投影する 物体の材質光の指定は，この例では object() の中でしてはいけない。両者とも同じ display() など表示 のレベルで行わないと，影の材質光が正しく出ない。影をぼかすには，他の複雑な方法を参考にしなければ ならないであろう。さらに，これまでの方法で得られる影は平面上の影であって，他のオブジェクトや自分自身への影は，オ ブジェクトを構成しているポリゴン上への影を計算し, 更に影のクリップの計算も生じる。
7
スプライン
7.1
はじめに
与えられた点を通過するまたは参照点とする曲線や曲面を描く場合に注意しなければならない問題を一 般論として考える．測量で得た有限個の点列を考えて，それらの点を順に通過する線路や道路を引く場合 を想定すればよい．最初から最後まで 1 つの関数で表すことはできないから，測定点を目標にして，順に 曲線で近似して継いで行くことになるが，継ぎ目で次のような点に留意する必要がある． 第 1 に，曲線や曲面は連続していなければならないことで，それぞれの端点での値が一致していなけ ればいけない．第 2 に，線路や曲線は滑らかに連続していなければならないことで，数学的には継ぎ目で 接線が一致するようにしなければならない．ただし，一致するといっても方向を含めて一致することを要 求する．即ち，接ベクトルが一致しなければならない．第 3 に，継ぎ目を通過する際に異常な力が掛から ないようにすることであり，力は加速度であるから 2 階の導関数が連続であることを要求する．例えば， f (x) =−x2 (x < 0), x2 (05 x)，という関数は x = 0 で両方の端点の値は 0 であるから連続してい る．次に，f0(x) =−2x (x < 0), 2x (0 < x) であるから x = 0 で共に 0 であり，2 番目の条件を満た している．しかし， f00(x) =−2 (x < 0), 2 (0 < x) であるから，x = 0 で 2 階微分の値がなく，3 番 目の条件を満たしていない．このような曲線の線路を走ると，継ぎ目の x = 0 を通過した途端に力の掛か る方向がいきなり逆方向になることを表しているから，この点で左右に激しく揺すられることになり，適切 な曲線でない．実際の線路や道路でも同じ配慮がなされる曲線を使用している． 以上のようなことに加えて CG ではレンダリングの速度を考慮しなければならないので，数学的に正確 であれば良いというだけでは不十分である．7.2
パラメトリック曲線
最初に多項式で曲線を描くことを考える．n + 1 個の点列 (x0, y0), (x1, y1), . . . , (xn, yn) を通る多項式 は n 次の多項式である．その係数の決定も連立一次方程式を解くだけで良いが，陽関数表示の性質から x0, x1, . . . , xn が単調列でない場合はこの方法は使えない．対策として，部分的に低次の多項式で表しなが ら継いでいくという方法が考えられるが，上記で説明したような接線を滑らかに継ぐという問題だけでも 解決できない．即ち，陽関数表示を使うことはこの場合不適当ということになる．次に考えられる方法は陰 関数表示 f (x, y) = 0 で曲線を表すことである．しかし，この表示法から曲線を描くには，高次の方程式を 解かねばならず，公式もアルゴリズムもない． したがって，数学では曲線とは 1 パラメータの関数の組 (f (t), g(t))，(f (t), g(t), h(t)) をいい，曲面とは 2パラメータ関数の組 (f (s, t), g(s, t), h(s, t)) をいうと定義している．CG ではこのようにパラメターで表 された曲線，曲面をパラメトリック曲線，パラメトリック曲面と呼んでいる．7.3
補間曲線
以下では多項式で表されるのパラメトリック曲線で考えよう．曲線を
(xi,0+ xi,1t +· · · + xi,mtm, yi,0+ yi,1t +· · · + yi,mtm, zi,0+ zi,1t +· · · + zi,mtm)
と表す．パラメトリック曲線の場合も，点列の初めから最後まで 1 つの式で表すのでなく，連続する点の間 を多項式で表し，式を変えながら点を継いでいくという考えを取る．全体を通して多項式の次数を同じに しておく． 始めに 1 次式で表されるか考えてみよう．1 次関数で表される曲線は直線であり，3 つ以上の任意の点を 通る直線はないので，折れ線近似しか得られない．折れ線は継ぎ目で微分不可能であるから 1 次関数は不 適当である． 次に 3 次式で表されるパラメトリック曲線を考える．どの座標成分でも同じ議論ができるので，代表し て 1 つの関数で考える．通過する点の列を (y0, y1, . . . , yn)とし，細分点列間 [yi, yi+1]で関数 fi(t) = ai+ bit + cit2+ dit3 (05 t 5 1) (11) を考える． 内部の点 y1, y2, . . . , yn−1 で考えると，上で説明した第 1 の条件を満たすために fi−1(1) = yi fi(0) = yi (12) 第 2, 第 3 の条件を満たすために fi0−1(1) = fi0(0) fi00−1(1) = fi00(0) (13) が i = 1, 2, . . . , n− 1 に対して成立する． また，端点 y0, ynで次の式が成り立つ． f0(0) = y0 fn−1(1) = yn (14) 以上で 4n 個の未知数に対し関係式が 4(n− 1) + 2 = 4n − 2 個であるから，残り 2 個の関係式があれば解 が一意的に求まる．次の 2 つの境界条件をつける．但し，一般には残り 2 つの条件は任意である． f000(0) = 0 fn00−1(1) = fn00(0) = 0 (15) 係数行列の独立性に関しては，異なる i の式と異なる階数の微分係数を使用した式であることより問題は 起きない．以上でこの連立方程式の解が一意的に求まる． 方程式は 1 次連立方程式ではあるが，未知数が 4n 個あり，区分曲線は互いに連結していることから，式 (12),(13)で見るように i 毎に分離できないので，非常に大きな行列の計算をしなければいけない．しかし， 次の簡潔で美しい計算方法が見つけられた．(Bartels et al., 1998) fi0(0)を使うと以下のように fi0(0)の漸化式が成立する．上の関係式を係数で表してみると，第 1 の条件 から fi(0) = ai= yi fi(1) = ai+ bi+ ci+ di= yi+1 (16)
第 2, 第 3 の条件から
fi0(0) = bi= Di fi0(1) = bi+ 2ci+ 3di= Di+1 (17)
この 4 つの関係式から次の漸化式式が得られる．
ci = 3(yi+1− yi)− 2Di− Di+1 (18)
di = 2(yi− yi+1) + Di+ Di+1 (19)
追加した境界条件 f000(0) = c0= 0 より 2D0+ D1= 3(y1− y0)を得る．次に f000(1) = f100(0) の条件より c1= 3d0，これに上で求めた式を代入して D0+ 4D1+ D2= 3(y2− y0)，以下同様の計算を続けて次の関 係式を得る． 2 1 1 4 1 1 4 1 .. . . .. . .. . .. 1 4 1 1 2 D0 D1 D2 .. . Dn−1 Dn = 3(y1− y0) 3(y2− y0) 3(y3− y1) .. . 3(yn− yn−2) 3(yn− yn−1) (20) 参照点が閉曲線の場合はこの式の最後の部分が次のようになる． Ã 1 4 1 1 1 4 ! Ã Dn−1 Dn ! = Ã 3(yn− yn−2) 3(y0− yn−1) ! (21) この行列の計算をすれば Di が求まり，Di を代入すれば ai, bi, ci, di が直ちに求まる． 以上のように未知数と方程式の個数の関係式から分かるように，2 次のパラメトリック曲線は決定不能で， 4次以上のパラメトリック曲線を決定するには，妥当な条件を多く付加しなけらばならないことが分かるで あろう．この意味で f000(0) = fn00(0) = 0 の条件をつけて得られた 3 次のパラメトリック曲線を 3 次スプラ
イン (cubic spline)，または自然なスプライン (natural spline) と呼ぶ．
