行列uniform(行列入門)

行列に慣れる

　慣れるために必要なことはとにかく触れることです。ピアノも毎日5分でいいから触るだけで違うそうです。なので行列に対する苦手意識も触れることで解消できます。そういうわけで行列に触れて苦手意識を解消しましょう。得意な人は、何でこんな簡単なことを...って思いながら読むことになると思います。それで構わないので、最後までお付き合いください。glslの行列の取り扱いに慣れるのが今回の主目的です。
　ていうか行列なんかどうでもいいからピアノ練習したい...

コード全文

作品14のリンク

プログラムの概要

　2次元のベクトルをひとつ用意します。2x2行列をひとつ用意します。ベクトルは列ベクトルです。行列の成分は行ベースで読みます。つまり、 \[ v=\begin{pmatrix} 0.1 \\[4pt] -0.1 \end{pmatrix},~~~~M=\begin{pmatrix} 2 & 3 \\[4pt] 2 & -1 \end{pmatrix} \] とあるわけですが、行列Mは2, 3, 2, -1と読みます。それで、これをvに作用させます。 \[ Mv = \begin{pmatrix} 2 & 3 \\[4pt] 2 & -1 \end{pmatrix} \begin{pmatrix} 0.1 \\[4pt] -0.1 \end{pmatrix} = \begin{pmatrix} -0.1 \\[4pt] 0.3 \end{pmatrix} \] 行列の掛け算やベクトルへの作用については特に説明しません。一応書いておきます。前提知識ですが、無いと困るので。 \[ \begin{pmatrix} a&b \\[4pt] c&d \end{pmatrix} \begin{pmatrix} x \\[4pt] y \end{pmatrix} = \begin{pmatrix} ax + by \\[4pt] cx + dy \end{pmatrix} \] \[ \begin{pmatrix} a&b \\[4pt] c&d \end{pmatrix} \begin{pmatrix} a'&b' \\[4pt] c'&d' \end{pmatrix} = \begin{pmatrix} aa'+ bc' & ab' + bd' \\[4pt] ca' + dc' & cb' + dd' \end{pmatrix} \] それで、これと同じことをプログラム内で実行します。得られる結果を、あらかじめ外で計算した結果と比較して、一致すればOK！というわけです。それだけのプログラムです。
　その計算方法をいくつか用意してコンフィグでいじれるようにしました。正解はすべて同じ結果で、正方形の位置と線をクロスさせる位置が一致します。不正解もすべて同じ結果で、この場合線をクロスさせる位置（＝正解の位置）と正方形の位置がずれます。つまり、不正解です。

uniformsの取得機構

　今回は正方形の描画と線の描画で2つのプログラムを使います。それゆえ今までのやり方では非常に不便です。なのでuniformsの取得処理をメソッド化しました。最終的にプログラムに紐付けます。 

// uniform生成関数
function registUniforms(gl, pg){
  // uniformの取得
  const numActiveUniforms = gl.getProgramParameter(pg, gl.ACTIVE_UNIFORMS);
  const uniforms = {};
  for(let i=0; i<numActiveUniforms; i++){
    const uniform = gl.getActiveUniform(pg, i);

    const location = gl.getUniformLocation(pg, uniform.name);
    uniform.location = location;
    uniforms[uniform.name] = uniform;
  }
  pg.uniforms = uniforms;
}

　今まで書いてきたコード、そっくりそのままです。違うのは最後にpgにuniformsをセットしているところだけです。これでプログラムだけでuniformの登録処理を実行できます。

行列uniform

　行列uniformの登録にはuniformMatrix[234]fvを使います。今回は2x2なので2fvしか使いません。レファレンスはこちらです。推奨されているのはFloat32Arrayですが、普通の配列でも問題ありません（列挙はもちろんダメ）。数でないものが入っていても特にエラーは出ませんがサイレントで不具合が起きます。まあFloat32Arrayが無難ですね。vと付いていますし配列指定なんですが、ユニフォーム名は配列でなければ普通に（この場合は）uMatです。配列の場合はuMat[0]などとなるのもこれまで通りです。

　（追記...通常の配列でもいいですが、その数はfloat,つまり4バイトの小数であることが想定されます。つまりどんな配列であってもFloat32Arrayとしてみなされるということです。なので小さい数であってもBYTEやSHORTのように扱うことはできませんし、場合によっては丸め込みが発生します。今後、たとえばindexBufferのときなんかにその違いを実感することになるかと思います。）

　それで、... 

  const mat = new Float32Array([2, 3, 2, -1]);
  const test = new Float32Array([0.1, -0.1]);
  const answer = new Float32Array([0, 0]);
  answer[0] = mat[0] * test[0] + mat[1] * test[1];
  answer[1] = mat[2] * test[0] + mat[3] * test[1];
  // 誤差は仕方ない。
  console.log(answer[0], answer[1]); // -0.1, 0.3

  gl.useProgram(pgCalc);
  gl.uniform2f(pgCalc.uniforms.uTest.location, test[0], test[1]);
  gl.uniformMatrix2fv(pgCalc.uniforms.uMat.location, false, mat);

　ここでは行列「mat」とベクトル「test」を用意しています。計算の結果は「answer」として計算で出しています。そのあと登録してますが、行列については2番目の引数をfalseにしていますね。レファレンスでtransposeとなっているものです。
　レファレンスにもありますがここは必ず「false」にしてください。レファレンスが「Must be（やるな！）」と言っているくらい絶対です。このプログラムではやりませんが、実はここをtrueにすると正解と不正解の内容がひっくり返ります。つまり不正解がすべて正解になり、その逆も然りです。それが何を意味するのかはこの後の説明を読めばわかると思います。ここではこれ以上触れません。
　そのあとの、これは正解の点を表す交叉線の描画用のuniformです： 

  gl.useProgram(pgAnswer);
  gl.uniform2fv(pgAnswer.uniforms["uAnswer[0]"].location, new Float32Array([answer[0], -1, answer[0], 1, -1, answer[1], 1, answer[1]]));

　uniformはプログラムを走らせないと登録できません。プログラム内で固定されるuniformは1回放り込めば十分なのでここで入れています。重複を防ぐのは負荷軽減というよりは単純に「不必要」だから、ですね。特に余計なことは考えていません。
　それでは、シェーダー内部での計算に移りましょう。

シェーダーでの行列計算

　描画関数はこちらです。黒でクリアして、このuPatternというのは計算方法の指定です。6パターン用意しました。0,1,2が正解で3,4,5が不正解です。結果となる位置に正方形が描画されます。正しければ、そのあとで描画する2本の直線がこれを貫く仕掛けです。 

  const drawResult = () => {
    gl.clearColor(0,0,0,1);
    gl.clear(gl.COLOR_BUFFER_BIT);
    gl.useProgram(pgCalc);
    gl.uniform1i(pgCalc.uniforms.uPattern.location, config.uPattern);
    gl.drawArrays(gl.TRIANGLE_STRIP, 0, 4);
    gl.useProgram(pgAnswer);
    gl.drawArrays(gl.LINES, 0, 4);
    gl.flush();
  }

正解パターン

　次の3パターンです。

  // 正解パターン
  if(uPattern == 0){ p += (v * uMat); }
  else if(uPattern == 1){ v *= uMat; p += v; }
  else if(uPattern == 2){ v = vec2(dot(v, uMat[0]), dot(v, uMat[1])); p += v; }

　まず1つ目ですが、vにuMatを右から掛けています。え？ですね。まあそれはそうです。CPUサイドでは左から掛けているので。しかし2つ目を見てください。vに右から「*=」でuMatを掛けています。とても自然です。glsl内部では基本的に作用素は右からやってくるので、外で列ベクトルに作用させるようにはどうしたって書けません。なので「右から」です。

　3つ目ですが、これは内積表現です。uMat[0]には初めの2つの成分が入っています。uMat[1]にはそのあとの2つの成分が入っています。実はこれらはvec2です。すなわち、uMatが[2,3,2,-1]で、uMat[0]がvec2(2,3)で、uMat[1]がvec2(2,-1)というわけです。なお、0番パターンを： 

    mat2 mm = mat2(uMat[0], uMat[1]);
    //p += (v * uMat);
    p += (v * mm);

こう書いてもいいし、

    mat2 mm = mat2(uMat[0].x, uMat[0].y, uMat[1].x, uMat[1].y);
    //p += (v * uMat);
    p += (v * mm);

こう書いてもいいです。つまり、外でも中でも行列の定義方法は基本的に同じです。まあ当然ですね。だって、外でも中でもやってることは、行列を用意して、ベクトルを用意して、それを行列に掛ける、ただそれだけですから。そして行ベース。行ベースです。
　とはいえ中でベクトルに行列を掛ける方向が逆なので、まあ解釈としては、こうですね： \[ \begin{pmatrix} 0.1 & -0.1 \end{pmatrix} \begin{pmatrix} 2 & 2 \\[4pt] 3 & -1 \end{pmatrix} = \begin{pmatrix} -0.1 & 0.3 \end{pmatrix} \] 全く同じ計算です。行列は転置されていますが、ベクトルに作用させるためにそう解釈しているだけです。ベクトルが横なのは、数学的に考えて2x2を右から掛けるなら縦だと具合が悪いからです。外でも、中でも、数学に従うのが常に自然な立場です。混乱を防ぐことができます。

　あとでちょっと述べますが、あまり転置に気を取られない方がいいです。「作用する」と考えてください。行列は線形変換です。線形変換Mとベクトルvがあって、Mがvに作用してベクトルv'が得られる、という「現象」そのものは中でも外でも一緒です。作用する方向が違うだけです。

間違いパターン

　ですから、次の3つは間違いということになります：

  // 間違いパターン
  else if(uPattern == 3){ p += (uMat * v); }
  else if(uPattern == 4){ vec2 w = uMat * v; p += w; }
  else if(uPattern == 5){ p += (v * transpose(uMat)); }

　演算上、実は行列をベクトルに掛ける処理は逆も可能です。数学的にはあり得ないですが、可能なことになっています。その意味ですが、パターン5を見てください。transposeを使うと転置を取ることができます。まあそういうことですね。これは転置行列を先ほどの内容で作用させる処理です。当然、全く異なる結果になります。作品サイトで確かめてみてください。
　一応書いておくと、 \[ \begin{pmatrix} 0.1 & -0.1 \end{pmatrix} \begin{pmatrix} 2 & 3 \\[4pt] 2 & -1 \end{pmatrix} = \begin{pmatrix} 0.0 & 0.4 \end{pmatrix} \] ですね。それっぽい位置になっているでしょうか。
　なお、掛ける方向をこの方向にする場合、当然ですがパターン1のような乗算代入記法は異なる結果を生じるため使えません。個人的には不便だと思います。まあ人の勝手ですが...

　単独は以上です。次に行列同士の掛け算が出てくるパターンをやりましょう。ここまでの議論が理解できるなら簡単です。

行列の掛け算の内部処理

　次に、複数の行列を用意して順繰りに掛けます。CPUサイドでベクトルに2つの行列を順繰りに（もちろん左掛けで）作用させて、内部でも同じことをして比較します。結合律があるので実際には行列の掛け算の結果を適用します。中では順繰りとかいろいろやります。なぜなら焦点は内部処理だからです。CPUサイドはどうなるかわかっているので、どう計算しようが何の問題もありません。
　じゃあさっそくやりましょうか。たのしみ！

コード全文

作品15のリンク

行列の配列uniform

　行列を3つ用意します。先ほどのやつと、新しいひとつ、さらにそれらの積です。これらをm0,m1とすると、先に作用するのはm0なので、形としてはm1*(m0*v)です。ゆえに(m1*m0)*vであり、vに作用するのは(m1*m0)です。自然な行列の積です。だからこれをm2として用意します。 

  const mat0 = [2,3,2,-1];
  const mat1 = [1,2,-2,-1];
  const mat2 = new Array(4);
  // mat1 * mat0 = mat2.
  mat2[0] = mat1[0]*mat0[0] + mat1[1]*mat0[2];
  mat2[1] = mat1[0]*mat0[1] + mat1[1]*mat0[3];
  mat2[2] = mat1[2]*mat0[0] + mat1[3]*mat0[2];
  mat2[3] = mat1[2]*mat0[1] + mat1[3]*mat0[3];
  console.log(mat2); // 6, 1, -6, -5

　一応検算します。 \[ \begin{pmatrix} 1&2 \\[4pt] -2&-1 \end{pmatrix} \begin{pmatrix} 2&3 \\[4pt] 2&-1 \end{pmatrix} = \begin{pmatrix} 6&1 \\[4pt] -6&-5 \end{pmatrix} \] 　ミスがあるといけないので。ミスしてると思ったらご報告ください。もちろん「m0*m1の方が正しいんじゃないか」とか意味不明なことを言われても困りますが...で、これらをまとめて登録します。 

  const mat = new Float32Array(mat0.concat(mat1).concat(mat2));
  const test = new Float32Array([0.1, -0.1]);
  const answer = new Float32Array([0, 0]);
  // もうダイレクトにmat2を掛けてしまう
  answer[0] = mat2[0] * test[0] + mat2[1] * test[1];
  answer[1] = mat2[2] * test[0] + mat2[3] * test[1];
  // 誤差は仕方ない。
  console.log(answer[0], answer[1]); // 0.5, -0.1

　CPUサイドなので素直に結合律に従ってm2を作用させてしまいましょう。行列は配列で登録します。行列は配列の場合も同じ関数を使います。配列で用意するのも一緒です。もちろんアクセスする名前は[0]が付きます。 

  gl.useProgram(pgCalc);
  gl.uniform2f(pgCalc.uniforms.uTest.location, test[0], test[1]);
  gl.uniformMatrix2fv(pgCalc.uniforms["uMat[0]"].location, false, mat);

　いいですね。今回は正解を6パターン用意しました。さっそく見ていきましょう。

シェーダー内の行列掛け算

　中でもやることは同じです。すなわち、ベクトルvにまずm0を掛けて、次にm1を掛けます。m0を掛けてから、m1を掛けます。それではまず正解のパターンを確認します。

正解パターン

  // 正解パターン
  if(uPattern == 0){
    p += (v * uMat[0]) * uMat[1];
  }else if(uPattern == 1){
    p += (v * (uMat[0] * uMat[1]));
  }else if(uPattern == 2){
    p += (v * uMat[0] * uMat[1]);
  }else if(uPattern == 3){
    v *= uMat[0]; v *= uMat[1]; p += v;
  }else if(uPattern == 4){
    v = vec2(dot(v, uMat[0][0]), dot(v, uMat[0][1]));
    v = vec2(dot(v, uMat[1][0]), dot(v, uMat[1][1]));
    p += v;
  }else if(uPattern == 5){
    p += (v * uMat[2]); // trueにすると9の結果になる
  }

　初めの4つのパターンはいずれも、m0を掛けてからm1を掛けるという処理を素直に表現したものになっています。3なんか一番露骨にそれをやっています。全部同じ結果です。掛ける方向が逆ですから、当然ですが、合成してできる行列は「"m0"*"m1"」という形になります...が、これは自然でしょう。作用する方向が逆である以上、この順でないと先にm0が作用できないからです。行列としては「m1*m0」です。おかしなことは何も起きていません。

　パターン4も基本的には同じことをしています...というかさっきのコードの処理を順繰りに実行しているだけですね。もちろん行列同士の掛け算をベクトル分解で書いてもいいんでしょうが、無駄に煩雑になるだけですし、どうしたって転置処理が絡むので、混乱を招くのでやめました。

　パターン5は配列に放り込んでおいた行列の積を使っています。もちろん同じ結果になります。先ほどのコードのようにバリエーションを用意しても良かったんですが、普通に冗長なので省略しました。

　なおパターン2ですが、まあこれでもいいんですが、可読性が落ちるのでお勧めしません。分けた方が分かりやすいと思います。好みの問題ですね。

間違いパターン

else if(uPattern == 6){
    p += (uMat[1] * (uMat[0] * v));
  }else if(uPattern == 7){
    p += ((uMat[1] * uMat[0]) * v);
  }else if(uPattern == 8){
    p += uMat[1] * uMat[0] * v;
  }else if(uPattern == 9){
    p += (uMat[2] * v); // trueにすると5の結果になる
  }

　外と同じ書き方をしました。まあそもそもこんな書き方ができることが問題なんですよね...数学的にあり得ないので。要するに転置とみなせば同じ結果になるんですが、面倒なので省略しました。

　ただ今回は行列の積が絡んでいるので若干異なる状況になっています。パターン9を見てください。これの結果は、まあ試せばわかるんですが、画面外に飛び出します。実は(1.2,0.6)が正解です。どういう計算かというと、こう： \[ \begin{pmatrix} 0.1 & -0.1 \end{pmatrix} \begin{pmatrix} 6 & 1 \\[4pt] -6 & -5 \end{pmatrix} = \begin{pmatrix} 1.2 & 0.6 \end{pmatrix} \] 行列の積は転置すると順序が逆になるため、全く違う行列が生成されます。それでこのような結果になるわけですね。
　コメントにも書きましたがfalseをtrueにするとこれの結果はパターン5になります。そしてパターン5の方がこれになり、吹っ飛んでいきます。6,7,8は0,1,2,3,4と入れ替わります。やらなくていいです。

転置と考えない方がいいよ

　転置と考えたい気持ちは分かります。mdnのサイトにもglsl内部では行列は列ベースで成分を書く的なことが一部にほのめかされています。しかしそれでもやはり数学的には行ベースで扱った方が自然ですし、列ベクトルに作用した方が自然です。そのような立場に固執する理由は、その方が、遍く存在する数学リソースを使いやすいからです。ほとんどの数学リソースはベクトルを列で扱っており、行列も左から掛けられています。それを再利用したい場合、列ベースで行列を考えるのは普通に不利です。頭が混乱します。特に顕著なのが、ベクトルの軸周りの回転です。

　もう一つの理由は、行列の転置は基本的には作用の結果だからですね。転置とは基本的に新しい行列を生み出すための処理です。実際、3Dの処理でとある過程で転置を取る必要が出てきます（法線など）。転置は基本的に新しい行列を生み出すための処理であり、記法上の都合による「見せかけの転置」とごっちゃにすべきではありません。あくまで作用する方向が逆であるということにのみ留意して、転置しているというイメージを極力持たないようにした方がいいと思います。

　行列は計算を鵜呑みにすると痛い目を見るので、思考停止する方向に誘導するリソースが多いのは分かるんですが、自分できちんと理解することをおすすめします...

行か、列か

　Qiitaの方ではもっと過激なことを書いてるんですが、こっちではほどほどに済ませて、次に行きたいと思います。
　行列の枠組みを作るうえで重要なのは「自分が混乱しないこと」です。行列は非常に強力です。行列が苦手な人には申し訳ないですが、3Dは行列が無いと何にもできないというのが厳然たる事実です。なのできちんと枠組みを作らなければなりません。その際に数学のマナーを守る姿勢は混乱しないようにするうえで非常に重宝します。なのでなるべく推奨したいわけです。
　それで、p5はこんな感じで書いてます。 

void main(void) {

  vec4 viewModelPosition = uModelViewMatrix * vec4(aPosition, 1.0);
  gl_Position = uProjectionMatrix * viewModelPosition;

  vec3 vertexNormal = normalize(uNormalMatrix * aNormal);
  vVertTexCoord = aTexCoord;
  ...
}

wgld.orgもこんな感じで、要するに本記事で述べたような書き方をしているのは少数派です。

attribute vec3 position;
uniform mat4 mvpMatrix;

void main(void){
    gl_Position = mvpMatrix * vec4(position, 1.0);
}

　頭が混乱するのでやりたくないんですが、当たり前ですが、次の計算とは異なる結果になるわけです： 

attribute vec3 position;
uniform mat4 mvpMatrix;

void main(void){
    vec4 p = vec4(position, 1.0);
    p *= mvpMatrix;
    gl_Position = p;
}

　理解が足りない場合、このように書き直すことで当然バグが発生するわけですが、その理由を理解することはおそらくできないでしょう。無駄な時間を過ごすことになります。

　WebGL Fundamentalsもこんな感じなんですね...もうこれがスタンダードになっているんでしょうね。 

attribute vec4 a_position;
 
uniform mat4 u_matrix;
 
void main() {
  // positionを行列に掛ける。
  gl_Position = u_matrix * a_position;
}

　左から掛けた方が分かりやすいのは分かるんですが、乗算代入と値がズレるのはどうしたってバグの原因になりやすいので推奨できないです。で、それによる齟齬が明らかに見て取れるのが、回転行列です。さっきのサイトで、2次元の回転では、 

  rotation: function rotation(angleInRadians) {
    var c = Math.cos(angleInRadians);
    var s = Math.sin(angleInRadians);
    return [
      c,-s, 0,
      s, c, 0,
      0, 0, 1
    ];
  },

となっており自然なんですが、3次元のz周りの回転、要するにこれに相当する回転ですが、 

  zRotation: function(angleInRadians) {
    var c = Math.cos(angleInRadians);
    var s = Math.sin(angleInRadians);
 
    return [
       c, s, 0, 0,
      -s, c, 0, 0,
       0, 0, 1, 0,
       0, 0, 0, 1,
    ];
  },

こうなっています。符号が、逆です。まあこれによりあのシェーダーでも正しく動作するわけですが、明らかに数学的不整合が起きています。この行列をCPUサイドで列ベクトルに作用させたら逆回転になってしまいます。CPUサイドでも行ベクトルで扱えばいいんですがそれだとリソースの活用が難しくなります。色々大変です。
　この理論には抜け道があります。「動けば何でもいいじゃないか」です。自分もそういう立場だったことがあるので大きなことは言えません。今は数学の方が大事なので、もうそういう立場は捨てました。理解が深まれば、最終的にはそっちの方が自由度は高くなります。結局、理解度がすべてです。

　なおUnityや他の環境では左掛けを推奨していたりする場合もあるそうですが、どのみち土台の（数学の）理解がきちんとしていないと仕様に振り回されて損することになるので、あやふやな理解だとおぼつかなくなるのはどこでコードを書くにしても同じかと思います。

　もちろんこれらは自分の利用する枠組みではないので、どんな仕様であろうとどうでもいいことです。大事なのは利用者自身が混乱しないことです。それさえ守られれば何の問題もありません。
　なお、この辺の理解は自分もずっと不十分で、fisceの前身のライブラリでは普通にベクトルに左から行列を掛けてました。記事を書くまで気にも留めませんでした。そういう感じなので、宗派替えした今では全く同意できないんですが、そういうやり方がしたい人がいることは自分もよくわかります。自分の好きな道を進んでください。

　行列については以上です。お疲れ様でした。

今回登場した関数

uniformMatrix[234]fv