2014年3月27日木曜日

フォトンマップの実装

イントロダクション

いままで本ブログでは、3DCGのレンダリング手法であるレイトレーシングとパストレーシングの記事をいくつか書いてきました。レイトレーシングはフォンシェーディングなどのレンダリング手法からすれば、反射・屈折などが扱え、格段に写実的になり高品質な画像を合成できるようになっています。ただし、レイトレーシングにも弱点がありました。フォンシェーディングはリアルタイムレンダリングに使える手法ですが、レイトレーシングでは、リアルタイムに描画するのは難しいでしょう。また、間接光や集光模様などの影響が扱えません。全体的に暗めの画像が生成され、影も違和感があります。写真のような画質とはとても言えません。そんなに時間をかけずに、そこそこの品質の画像を生成するのにはよいでしょう。写真のような画質の画像を生成するには、すべての光の軌道を計算する必要があります。現実的にその計算を行うことは不可能ですが、それに近いことを行うのがパストレーシングです。パストレーシングは、全ての光の軌道を計算する代わりに、その中のいくつかの光の軌道をランダムに選び、全体像を推定する統計的な手法です。選ぶ光の軌道の数が少なければ、正確に推定することができずに、生成された画像にムラができてしまいます。なので、十分な数だけ光の軌道を選ぶ必要があります。こうして生成された画像は、写真とは見分けがつかないほどリアルなものになります。ただ、十分な数の光の軌道を追跡するのは、かなり時間がかかるのが欠点と言えるでしょう。本記事では、それらの欠点を補うために、フォトンマップを使った手法を取り上げます。


図1. レイトレーシングとパストレーシングの比較

パストレーシング

フォトンマップに行く前に、パストレーシングの問題点を振り返っておきましょう。パストレーシングでは、まず視点(カメラ)からスクリーン上のあるピクセルに向かって光線を放ちます。その光線はシーン中の視点にもっとも近いポリゴン上のある1点と衝突します。その表面の材質によって光線のその後の軌道が変化します。鏡面や屈折面だった場合は、正確にその後の軌道が計算されます(フレネルの法則)。拡散面(ざらざらした面)だった場合は、光は乱反射します。乱反射するということは、その後の光の軌道は、様々な方向へ光は分散していくということです。


図2. スクリーンのあるピクセルに放った光線


図3. 鏡面反射と乱反射

ちなみに、分散レイトレーシングでは、拡散面に衝突した時は、実際に複数(例えば、100本だったり1000本だったり)の光線を追跡します。ただ、この方法だと拡散面に衝突する度に、100本や1000本の光線が生まれ、結果的に膨大な光線(100のn乗や1000のn乗、nは拡散面との衝突回数)を追跡する羽目になります。

一方、パストレーシングでは、光線が拡散面に衝突した時、考えうる反射方向の中から、1本だけランダムに軌道を選びます。したがって、視点から放たれた光線は、最初から最後まで1本です。拡散面では、光線の反射方向はランダムに選ばれるので、状況によって光線の反射方向は毎回違います。つまり、スクリーン上のあるピクセルの色は、視点から放たれた光線の経路によって変わってくるので、もしかしたら、ある時は赤、またある時は青といったように状況によって(乱数生成環境によって)変わってきます。

この状況を現実世界に置き換えて考えてみると、光源から放たれた1つの光(光子)が、いろんな経路を通って目に届くことを意味しています。パストレーシングは(レイトレーシングも)、光源から来た目に届くはずの光を、逆に目から光源に向かって追跡しているのです。別々の経路を通った光は、様々な材質の表面を通ってくるので、やはり、それぞれ色(波長)が違うものになるはずです。ただ、現実世界では、光源から放たれた光線は、たったの1本だけが目に届くわけではなく、複数の光線(それこそ無限近いの光線)がいろんな経路を通って目に届くはずです。目に届く光は、無数の光線の合計の色として認識されます。

パストレーシングも現実世界の例に習い、光線の方向は逆方向ではありますが、視点から1つのピクセルに向かって複数の光線を放ちます。ある光線は、いくつかの反射・屈折の後、光源まで届くでしょう。また、ある光線は、十分な回数の反射・屈折を経た後でも、光源に到達することはないでしょう。実際に、光源から目に届く複数の光の経路というのは、視点から放った光線の内、光源に到達した光線の経路と同じになるはずです。つまり、パストレーシングは、十分な数の光線を使用すれば、ほぼ完全に現実世界をシミュレートすることができます。

しかし、十分な数の光線を使用しなかった場合は、合成された画像には、ノイズが生じます。このノイズを説明するために、サイコロを振ったときに出る目を例にあげましょう。例えば、サイコロを5回振ったとき、出た目が順に、3、2、2、1、2だったとしましょう。この場合、出た目の平均は、(3 + 2 + 2 + 1 + 2)/5 = 2 になります。もし十分な回数だけサイコロを振っていれば、おそらく平均値は3.5になるはずです。ある程度の回数まで試行すれば、平均値が3.5からあまりズレなくなるでしょう。ですが、試行回数が少ないと平均値が3.5から大きくずれることが多いでしょう。 パストレーシングにおいても、同様のことが言えます。視点からあるピクセルに放出する光線の数が少ないと、実際の色から大きく外れた色になる確率が高くなります。光線の数は十分に大きくとらないといけないのです。


図4. ノイズの生じたパストレーシング

パストレーシングの最大の難点は、この部分にあると言えるでしょう。十分な数の光線を使用すると、描画に膨大な時間がかかります。しかも、いくつかの場合は、ほとんどの光線が光源に到達することがありません。光源に到達しない光線というのは、無駄に時間を消費してしまっており、パストレーシングによるレンダリングの時間がかかる原因の1つと言えるでしょう。もっと効率的に、計算する必要がありそうです。

この問題を解決しようと試みるのが、本記事の主題であるフォトンマップを使用したレンダリング手法です。この手法を使えば、パストレーシングに比べ画像生成は格段に早くなります。ただし、近似計算的要素が含まれているため、完全な現実世界の再現にはなっていない部分もあると思いますが、その影響は小さくすることができます。

フォトンマップ

パストレーシングの問題点を解決するには、拡散面で反射した光線が光源に到達するまでの経路探索を、より効率よくすることです。そこで、拡散面上のある点Pにおける光の入出力に注目して、人の目に届く光の量を考えてみましょう。

点Pを経由して目に届く光の量は、光源から点Pに届く光の量に依存しているはずです。ここで、点Pに関する光の量を2つに分類します。1つめは、光源から点Pに届く光の量。2つめは、点Pから目に届く光の量。まず、2番目に分類した光の量について考えましょう。つまり、点Pに届いた光の量のうち、どの程度が目に届いているのかを考えます。完全な拡散面の場合は、点Pに入射した光の方向に依らずに、どの方向にも一定の割合で反射するので、点Pからは一部の光だけが目に届きます。また、目に届く光の量は、点Pから目までの距離に依存します。例えば、点Pに入射した光子の個数が100個だったとすると、反射して出て行く光子の数も100個です。そして、この100個の光子は、この後どこにいくのか? 簡単の為に、この光子は1秒間の間に1m進むものと仮定すると(本当は秒速30万キロ)、1秒後には、100個の光子は、点Pから1mだけ離れた場所にいるはずです。ある点から定距離だけ離れた点の集合が作る曲面は、球面です。なので、100個の光子は1秒後には、半径1mの半球面上(半球なのは、光子は裏側へは反射しないから)に一様に分布しているはずです。この時の球面の面積は、球面積の公式を使って、4 x π x 1 x 1 = 4πの半分の2πです。光子の数の密度は、100/2πとなります。点Pの光の密度の1/2πになります。光の明るさは、目に入ってくる光子の量に比例するので、点Pのすぐ近くにいる人と点Pから1mだけ離れている点にいる人が感じる光の明るさの比は、2π:1となります。では、rメートルだけ離れた人が感じる明るさは、どの程度でしょうか? 点Pから射出された100個の光は、rメートル進んだ時には、半径rメートルの半球面上に一様に分布しているので、rメートルだけ離れた場所の光の密度は、100/4πr^2になります。これが光の明るさです。つまり、光の明るさは、距離の二乗に反比例して弱くなっていきます。余談ですが、重力や電磁気力(そもそも光は電磁場の振動)などの力にも同じことが当てはまります。これらの力も同じような理由(この場合はフォトンの密度ではなくて、球面状の場の密度)で距離の二乗に反比例します。

次に、光源から複数回の反射・屈折を経て点Pに入射してくる光の量を考えてみましょう。パストレーシングでは、この量を点Pから光源に向かって追跡して求めました。前述したように、点Pからどの方向へ光線を飛ばせば、複数回の反射・屈折を経て、光源に衝突するのか分からないので、ランダムに光線を飛ばします。どこにも反射せずに、直接光源に届く経路を見つけるならまだしも、複数回反射・屈折して光源に到達する経路を見つけるのは、ランダムな方向に光線を射出していたのでは、偶然に衝突するまで待たなければなりません。あまり効率が良くないのです。しかも、現実に近い輝度を推定するには、十分な数の光線が光源に衝突しなければなりません。その十分な数の光線を得るには、膨大な時間がかかります。この点が、パストレーシングの難点でした。

パストレーシングのような方法とは別に、点Pに入射してくる光の量を推定する良い方法はないでしょうか。この問いに対する一つの解決案は、とても単純です。点Pから光源に至る経路を探索するのではなく、現実世界と同じように、光源から光線(光子=フォトン)を追跡すればよいのです。光源からフォトンを追跡すれば、全てのフォトンが輝度をもっているので、シーンのどこかで衝突すれば、それは意味のあるものになります。パストレーシングでは、放った光線が光源に衝突するのを待つのに比べたら、だんぜん効率がよさそうです。また、シーン全体にどのようにフォトンが分布しているかを知ることができるので、一度フォトンの分布を調べておけば、点Pの輝度を推定するだけでなく、他の点の輝度も推定できることがパストレーシングの場合と違います。パストレーシングでは、物体上の各点から、ランダムに光線を放って、光源までの経路を推定しなければなりませんでした。いろいろな場所で同じことを繰り返さなければならないのです。逆に光源から追跡する方法だと、一度だけシーン全体のフォトンの分布を計算すれば、輝度を推定したい他の場所でも使い回すことができます。この点は、大幅な計算量の削減に貢献するでしょう。このようなシーン全体のフォトンの分布はフォトンマップと呼ばれています。では、なぜ今までこの安直な方法(誰でもまず最初に考える)を試してこなかったのでしょうか?それは、光源からフォトンをランダムに放っても、点Pに届く確率がとても小さくて、非常に効率が悪いからだと思います。点Pには、ほとんど面積がないので、フォトンが点Pとほとんど衝突しません。これでは、点Pに入射する光量を推定することが困難です。ただ、これが点Pの周辺まで含めれば話は変わります。点Pの周辺も点Pとだいたい似たような輝度だろうという大胆な仮定を取り入れるならば、点Pの周辺(点Pを中心とする球)に入射する複数のフォトンから点Pの輝度を計算することができます。ただし、点Pの周辺は大体同じような輝度だという仮定は、常に成り立つわけではありません。壁と床の境目などの場所では急激に輝度が変化するので、そのような場所にはこの方法は、うまく働きません。


図5. 光源からフォトンを放つ・点Pの周辺のフォトンを収集

とりあえず、ここまでの方法をまとめましょう。
1. 光源からフォトンを放ち、シーン全体のフォトンの分布フォトンマップを構築する。
2. 視点からスクリーン上のあるピクセルに向かって放った光線と拡散面との交点をPとする。
3. 点P周辺のフォトンを収集し輝度を推定する。
4. 点Pの輝度から目に届く輝度を計算する。(2πr^2で割る)

フォトンマップの構築

物体上のある点(とその近傍)に入射する光の量を求めるには、フォトンマップを構築しなければなりません。フォトンマップは、シーンの中の拡散面に衝突する度にフォトンの情報を保持しています。保持する情報としては、衝突した位置、フォトンの色などです。フォトンの発生源はシーンにある全ての光源です。光源の種類によって違いますが、基本的に物理的に可能な方向へランダムに発射します。たとえば、平面光源の場合は、光源上からランダムに1つの点を選び、その点を中心に、天頂角が180度以内の方向に(半球を形作るように)ランダムにフォトンを発射します。一度拡散面に衝突したフォトンは、フォトンマップに格納され、そこで追跡をやめるのではなく、その後もフォトンを追跡していきます。

フォトンマップのデータ構造にも、注意を払う必要があります。というのも、ある点に入射する光の量を計算するときに、その点の近傍にあるフォトンを、フォトンマップの中から収集しなければならないからです。ただの配列にフォトン情報を詰めてもいいですが、ある点の近傍にあるフォトンを見つけるには、配列にある全てのフォトンとの距離を計算する羽目になります。探索が容易なデータ構造を選ぶのよいでしょう。ちなみに、今回、僕は八分木を使いました。


図7. 光源からシーンに放たれたフォトンの分布

レンダリング

フォトンマップを構築したら、レイトレーシングでフォトンを収集していきます。視点からスクリーン上のあるピクセルに向かって放たれた光線は、拡散面に衝突するまで反射屈折を繰り返します。拡散面に到達したら、その点の周辺のフォトンをフォトンマップから収集して輝度を計算します。点の周辺が意味するところは、さまざまあると思います。基本的には、その点を中心とした球の内側に入るフォトンを輝度の計算に使います。例えば、ある点中心にして、フォトンが100個収まるような球の使うという方法があります。100個が収まるには、ある点では、大きな球が必要かもしれませんし、別の点では、小さな球で済むかもしれません。いずれにせよ、各点で球の半径が変わります。そして、100個のフォトンが半球面を通過して、点Pの周辺に到達したことになります。基本的に、その点の入射輝度は、単位半球面を通過するフォトンの合計です。なので輝度を計算するには、球の大きさが単位半球だったら、フォトンがいくつ入射するかを推定しなければなりません。それは100個入射したときの球の表面積と単位球の表面積の比から求めることができます。例えば、100個のフォトンが収まる球の半径が10だったとしたら、表面積は4 x π x 10 x 10 = 400πで、単位球の表面積は、4πです。したがって、単位球の表面積は、半径が10の球の表面積の1/100になります。おそらく、入射するフォトンの個数も1/100になり、同時に輝度も1/100になります。なので、その点の入射輝度は、100個のフォトンの合計輝度の1/100です。これが各点における入射輝度の計算方法です。

別の入射輝度の計算方法として、各点で球の特定の数のフォトンの個数を収集できる大きさの球を使用するのではなくて、半径を固定して、その球の中に入ってくるフォトンで輝度を計算することもできます。どの点でも、同じ大きさの球を使用しているので、球に入って来たフォトンの合計の輝度の比がそのまま、各点の輝度の比と等しくなります。生成された画像に多少ムラができやすかもしれませんが、こちらの方法の方が簡単です。

先ほども述べましたが、壁と床の境なのどが輝度計算に使用する球の中に含まれてしまうと、正確に輝度を計算することができません。床のある点の輝度を計算しているのに、壁に到達したフォトンも輝度計算に入れちゃうからです。この現象を避けるには、ある点の近傍の範囲を決めるのを球を使わずに、円盤上のものを使うのがよいでしょう。輝度を計算したい点があるポリゴンの法線方向に球を潰したような立体です。こうすることで、法線と直交する平面(ポリゴン)に到達したフォトンを収集できる割合が高くなり、関係ない物体に届いたフォトンの収集を少なくすることができます。

フォトンが楕円体の中にあるか外にあるかの判定は、楕円体の中心とフォトンまでの距離と楕円体の中心からフォトンがある方向の楕円の表面までの距離を比較することで判定します。中心から表面までの距離は、
a * b / sqrt( b * cos(θ) * b * cos(θ) + a * sin(θ) * a * sin(θ) )
で求めることができます。aは長軸の長さ、bは短軸の長さ、θは法線と中心からフォトンまでの線分のなす角です。


図8. 円盤状の球体でフォトンを収集する


図9. 楕円を使わなかった場合、面と面との境付近で余分なフォトンが含まれている

このようにして収集したフォトンから、入射輝度(色)計算します。1つのフォトンが運ぶ光の量は同じなので、入射輝度は、単純にフォトンを足し合わせることで計算できます。入射輝度が計算できたら、その点から視点に届く光の量を計算します。この点は拡散面上の点なので、届いた光は一様に、様々な方向に反射して生きます。その一部が目に届くのです。この場合の目に届く光の強さは、先ほど述べたように1/2πr^2です。

あとは、スクリーン上の全てのピクセルに大して同じことを繰り返すことで、画像を生成することができます。

レンダリング結果

フォトンマップを使った方法では、光源から放射するフォトンの量によって、生成される画像の精度が変化します。その変化の様子の例を以下に示します。


100個フォトンを使用し、半径が0.5の場合。


1000個フォトンを使用し、半径が0.3の場合。


10000個フォトンを使用し、半径が0.1の場合。


100000個フォトンを使用し、半径が0.3の場合。


100000個フォトンを使用し、半径が0.1の場合。


500000個フォトンを使用し、半径が0.1の場合。


1000000個フォトンを使用し、半径が0.2の場合。

2014年3月17日月曜日

モンテカルロ積分を使ったレイトレーシング

本記事では、表面積分を統計的に行って光源からの寄与を計算するレイトレーシングを行います。積分を統計的処理で計算する手法をモンテカルロ積分というらしいです(本記事がそれに該当するかは不明。ただ、似たようなことをやっていると思います)。解析的に(数値計算や近似を使わずに式変形だけで正確に)積分計算をするのが難しい積分値を計算する時の手法です。ランダムに値を生成することから、カジノで有名なモンテカルロの名前を与えられたと聞いています。前回の記事では、スクリーン上のあるピクセルの色を計算するとき、光源からの寄与は、ある場所(輝度を計算したい場所)から光源を見たときの視野角に比例すると仮定して計算しました。今回は、もっと正確に計算するために、ある場所を中心とする単位半球面上に投影される光源の面積を統計的に求めます。つまり、光源がありそうな方向一帯に光線をいくつもランダムに放ち、光線が光源と衝突した回数と衝突しなかった回数との比から面積を求めます。もちろん、放つ光線の数が多いほど正確に面積を求めることができます。

面積の統計的推定

輝度を計算したい場所(視点から放たれた光線と視点に最も近いポリゴンとの交点)から見て、光源がありそうな方向の決定には、前回計算した視野角を使用します。ここで言う視野角は、図のように、視点から光源の外接円(または外接球)が収まる角度です。つまり、ある場所と光源の中心とを結ぶ線分と単位半球面との交点を中心として、半径が視野角の半分となるような単位半球面上に描かれる円の内側のみに光線を放ちます。円の内側のみ放つというのは、この円の面積は解析的に求めることができるので、光源がないとわかっているところに光線を放っても意味ないからです。今回は100本の光線を放ち面積を推定しました。100本のうち、ある光線は、光源にぶつかり、ある光線は、光線にぶつからないでしょう。そして、何本が光源に衝突するかがわかると、衝突した回数を100(放った光線の本数)で割れば図の単位半球面上の円の中に占める光源の割合が分かります。球面上の円の面積は、解析的に求めることができるので、その面積と先ほど求めた割合をかければ、求める光源の面積が計算できます。もちろん、光線の本数が多ければ多いほど正確に面積が求めることができるわけです。100本とか200本とか1000本とかで割合を計算してみて、結果にばらつきがないくらいの本数だったら十分だと思います。単位球面上の円の面積は、sinθdθdφを0≤θ≤'視野角/2', 0≤φ≤πの範囲で積分した値になります。
ここまでの過程で、ある場所の光源からの輝度が計算できます。輝度は、先ほど求めた光源の面積に比例するはずです。あとは、視点からスクリーン上の全てのピクセルに向かって光線を放ち、視点に最も近いポリゴン(ポリゴンでなくともよいが)との交点について、それぞれ同じ計算を施せば、1つの画像が生成されることになります。当然、前回に比べて計算量が増えるので、描画処理は遅くなります。

結果比較

前回の結果と比較してみます。前回は、単位半球面上に投影される光源の面積は視野角に比例すると仮定して輝度を計算しました。以下の最初の図(左側の図)が、今回行った統計的に面積を求めた時の結果です。次の図(右側の図)が前回の結果です。今回の結果の方が壁面の影のエッジが丸っこくなっている感じがします。おそらく、今回の結果の方が正しい結果なはずなんだけど、果たして実際にそうなんだろうか。現実は、間接光の影響があるので、直接光だけの影響を受けた物体の写真あるはずもないので確かめられないのだが。うーむ。

2014年3月16日日曜日

CygwinターミナルでCTRL+RIGHT/LEFT

Linuxのターミナルで単語単位でカーソルを移動したいとき、Ctrl+→Ctrl+←などのキーボードショートカットを使うことがあると思います。ただ、Cygwinターミナルのデフォルトの設定では、Ctrl+→Ctrl+←で単語単位で移動することができません。そんなときは、~/.inputrcに以下のような記述を追加すれば、Cygwinターミナルで単語単位の移動がCtrl+→Ctrl+←でできるようになるようです。

~/.inputrc
"\e[1;5C": forward-word   # ctrl + right
"\e[1;5D": backward-word  # ctrl + left

2014年3月1日土曜日

レイトレーシングの拡散面の輝度計算

3次元グラフィックの描画方法の一つであるレイトレーシング。光の反射・屈折などが正確に計算できる手法で結構キレイに画像生成ができます。ただし、間接光が計算できないので、暗めの画像になってしまいます。本ブログでも、Scalaで実装したレイトレーシングで描画した画像を載せていましたが、実は今までレイトレーシングよって生成される画像の各ピクセルの色の計算を結構適当に計算していました。今回は、この色の計算について考えてみたいと思います。

レイトレーシングの復習

簡単にレイトレーシングの仕組みについて復習してみましょう。現実の世界で、光源から発せられた光が目に届くまでには、いろんな物体に何回も衝突し、その一部が目に届きます。写真のような画像を作りたければ、光源が発せられた全ての光の道筋を計算すればよいことになります。でも、光源からの光の道筋を全て計算しても、最終的に目に届く光というのは、ごく一部です。ほとんど無駄な計算になります。コンピュータの使えるリソースは、限られているのに無駄な計算を何日もかけて計算する気にはなりません。もっと効率的に計算できる方法はないか。ちょっとぐらい写真画質じゃなくてもいいよ。と思った人が、たぶんレイトレーシングという効率的な方法を考えだしたんでしょう。

光源から光を追跡したら、膨大な計算量を必要とする上に、目に届くのは一部。じゃあ、目に届いた光だけを逆に追跡すればいいんじゃね。そう考えたんですね。目に届いた光というのは、スクリーン上のある点と目のある点を結ぶ半直線です。つまり、光源からの光を追跡するんじゃなくて、目がある位置からスクリーン上のある点を通過する半直線(光線=レイ)を追跡(トレース)しようってわけです。これがレイトレーシング(光線追跡)です。光源からの光は無数にありますが、目からの光線は、目がある場所からスクリーン上の各ピクセルとを結ぶ半直線のみなので、ピクセル数の光線を追跡するだけで済みます。

これで、目に届く光だけを計算することができます。目から出た光線は、シーンの中にあるいろんな物体にぶつかって、材質によって反射したり屈折したりして、ある光線は、光源にぶつかり、ある光線は、どこにもぶつからずに闇に消えていきます。光源にぶつかった光線は、その光線が通過したスクリーン上の点を明るくするでしょう。どこにもぶつからなかった光線は、その光線が通過したスクリーン上の点を明るくすることはありません。このようにして、各ピクセルの色が決められていきます。なので、物体表面の反射や屈折を計算して、どの光線が光源とぶつかるかを計算するのが重要になってきます。光線の反射や屈折の計算は、鏡や透明なガラスなどの材質だったら簡単です。では、ざらざらな表面の場合は、どうでしょう?ざらざらした表面は、光が乱反射します。ちなみに、このような性質の表面を拡散面と呼ぶようです。拡散面では、光線は、いろんな方向に反射します。鏡みたいに一方向に反射しません。ん?ちょっと待てよ。色んな方向に反射する?ということは、無数に分散する光線を追跡していかなきゃならんのですね。これじゃ、光源から発せられた光線を追跡するのと変わらんやん!だめやん。振り出しに戻ったやん。と思ってしまうかもしれませんが、レイトレーシングでは、拡散面にぶつかったら、そこで光線の追跡終えます。これがレイトレーシングが正確にピクセルの色を計算できないとこです。では、どうするのかというと、拡散面にぶつかった光線のピクセルの色は、拡散面の位置と光源の位置から近似計算します。間接光の計算はしません。間接光を計算するには、拡散面で分散した光線を追跡しなければならないからです。ちなみに、この追跡を統計的に行うのがパストレーシングです。ということで、レイトレーシングは、拡散面は近似計算で鏡やガラスなどの材質は、正確に計算するレンダリング手法です。

拡散面に届く光の量と目に届く光の量

先程も、拡散面では近似計算して、ピクセルの色を計算すると言いました。ここでは、具体的にその計算方法を考えてみます。以下のピクセル値の計算は、僕が勝手に近似計算しているだけなので、よりよい方法があるかもしれませんし、間違っているかもしれません。あしからず。

拡散面から来る光の量を計算を簡単にするために、光の量の種類を2つに分類します。1つは、光源から拡散面に届く光の量で、2つめは、拡散面から目に届く光の量です。

拡散面に届く光の量の近似計算

光源から拡散面上のある点に入射してくる光の量を正確に計算するには、その点を中心にした単位半球面を考えます。半球面は、その点での法線方向側にあります。つまり、拡散面の表に半球面があります。その単位半球面を横切る光線の数が、この点に入射してくる光の量です。これは、入射光線ベクトルを半球面に沿って積分することと同じことです。半球面に入射してくる光線ベクトルは、光源から直接入射してくるのもあれば、何回も反射してきて、入射してくるのもあります。前者の光線は直接光、後者の光線は間接光です。先程も言及したとおり、この入射ベクトルを全て計算するのは、大変なのでレイトレーシングでは、この積分計算を正確に計算せずに近似計算します。
まず、間接光は計算しません。間接光を計算するには、何回も反射した光線を考慮しなければならないので、非常に計算量を要するからです。直接光だけを考えましょう。半球面上を横切る直接光の数は、半球面上に投影される光源の面積に比例するはずです。半球面に投影される光源の面積は、図にあるように、魚眼レンズで見た時の光源が形作る図形の面積のようになると考えるとイメージしやすいでしょう。それでは、この面積を計算しましょう、と言いたいところですが、この積分計算もややこしそうなので、思い切って端折って、この直接校の計算も近似することにします。
では、どのようにして近似しましょう。たぶん、半球面上に投影される光源の面積は、図のように、拡散面上の点から光源を見た時の角度にだいたい比例するんだと思います。きっと。

以下にその角度を計算するコードを示します。
  val v0 = p -> center // 拡散面から光源の中心までのベクトル
  val d = radius // 光源を囲む円の半径
  val nt = normal * (v0 * normal) // v0の光源の法線方向成分
  val vr = v0 + (nt -> v0).normalize * d // 拡散面から光源の端までのベクトル
  val dt = vr angle v0 // 拡散面から光源を見た時の角度の半分

目に届く光の量の近似計算

次に、拡散面から目に届く光の量を考えます。拡散面に入射した光の量の一部が目に届きます。拡散面なので、この面に入射した光が全て目に届くはずがありません。そのごく一部の光が目にとどきます。正確には目に届く光は、ちょうど目のある位置の方向に偶然反射した光だけです。したがって、拡散面で反射する光が、どの方向にも一様に反射すると仮定すると、目に届く光の量は、拡散面に入射した光の量を半球面の面積で割った量だけになるはずです。ここで半球面の半径は、拡散面上の点から目の位置までの距離です。

それに対して、鏡で反射して目に入ってくる光は、鏡のその点に入射した光が全て同じ方向(目がある方向)に反射した光です。なので、この場合は、半球面の面積で割る必要がありません。

描画結果

レイトレーシングで描画した結果です。この画像をみる限り、光の量の計算は、よい近似になってるのではないでしょうか。