レイトレーシング(4): フォトンの生成
光源の定義をやりなおす
前回の記事で、光源について型クラスLightを定義した上で
各種光源の型(例えば点光源PointLight)をそのインスタンスと
するようにした。その後、シーン情報内で光源のリストを作ろうと
した時に問題に気づいた。Haskellでは"同じ型"しかリストに
できない。新たに面光源を用意したとして、「光源リスト」に
点光源と面光源を混ぜることはできない。
オブジェクト指向でのスーパークラスの概念とごちゃまぜに
してしまったのだ。。。
改訂した定義は以下の通り。
data Light = PointLight Color Flux Position3 flux :: Light -> Flux flux (PointLight _ f _) = f generatePhoton :: Light -> IO Photon generatePhoton (PointLight c _ p) = do : (あとは同じ)
引数でPointLightをパターンマッチさせるところがHaskellらしい。
今後、平面光源や球光源を追加するときは、dataの定義や各関数の
パターンを増やしていけばよい。
フォトンの生成
前回の話に立ち返り、メイン処理の前半を詳細化していこう。メイン処理を次に 再掲する。
main = do photons <- generatePhotons nphoton lgts photoncaches <- tracePhotons objs photons a <- forM (photoncaches) $ \i -> do putStrLn $ show i return ()
このphotons <- generatePhotons nphoton lgtsの部分だ。
nphotonとlgtsはそれぞれ、追跡したいフォトン数と光源のリスト。
次のように定義しよう。
nphoton = 100000 :: Int pl1 = PointLight (initColor 1 1 1) 2.0 (initPos 0 3.95 2.5) -- 2W lgts = [pl1]
フォトン数は生成される画像の品質を左右するパラメータだ。とりあえず10万個としよう。
点光源1個で試してみる。この点光源の定義は、赤緑青が同じ比率＝白色光
(initColor 1 1 1のところ)であり、発光の強さ(光束、ざっくり言えばエネルギー）が
2[W]、光源の位置は x,y,z=0, 3.95, 2.5 [m]ということ。ちなみにSI単位系とし、
長さは[m]とする。強さが2[W]だとかなり暗いのではと思われがちだが・・・。
光の強さというか明るさについては、白熱電球だと100[W]とか60[W]とかがメジャーだ。
ちなみに白熱電球の場合、エネルギー変換効率は2.3〜2.6%ぐらいだそうで、ほとんどが
熱などになってしまうということだ。地球に優しくないとして製造中止になるわけだ。
一方、最近LED電球ではlumen(lm)を明るさを示すために使っているが
白熱電球の100[W]クラスは大体1520[lm]ぐらいらしい。それで、lumenをエネルギーという
意味でワットに換算するとき、波長555[nm]の光では、
1[W] = 683[lm]
だそうなので、波長の違いをとりあえず無視して、1520 ÷ 683 ≒ 2[W] としたのだ。 この換算が正しいかどうかよくわからないが。 ただし、波長が違うとこの数字が大きく異なるので、本当は無視するとダメかもしれない。
ではあらためてgeneratePhotonsの定義を考えよう。要は、複数の光源からn個の
フォトンが放出される、という状況を作り出したい。Lightの定義で、ある光源から1個の
フォトンを放出する関数は定義した(generatePhoton :: Light -> IO Photon)。
これを使って全光源から放出される数をnphotonにしたいのだ。そうすると次のことを
考えていくことになろう。
ということで、コードにしてみた。
generatePhotons :: Int -> [Light] -> IO (Double, [Photon]) generatePhotons nphoton lights = do let power = (sum $ map flux lights) / (fromIntegral nphoton) ns = map (calcN power) lights photons <- mapM (mapM generatePhoton) (zipWith replicate ns lights) return $ (power, concat photons) calcN :: Double -> Light -> Int calcN power light = round (flux light / power)
ここで、powerがフォトン一個のエネルギー(W)でありnsは各光源が放出すべき
フォトン数のリストだ。
ここまでくれば後段はその数だけフォトンの生成を繰り返せば良い。具体的には各光源を
放出するフォトンの数だけ並べたリストを作り、そのリストの要素毎にgeneratePhotonを
呼び出しているだけだ。これで欲しい数のフォトンが得られた。
生成されたフォトンの品質
さて、これで10万個のフォトンを生成したが、果たしてちゃんとランダムに生成できている のだろうか?放出される偏りはないのだろうか? 点光源なので、無作為に球状のどの方向にも偏りなく放出されて欲しい。もし完璧に偏り なければ、フォトン数を無限に近づけることで全フォトンの放出方向のベクトル（大きさは 1に正規化されている）を足し合わせれば「0ベクトル」になるだろう。既出のメインルーチンの 処理を「途中下車」して、生成されたフォトン群を出力させた。さらにそれを下記の プログラムで読み込み、全ベクトルを足し合わせたものを出力してみた。
main :: IO () main = do np <- getLine pw <- getLine let nphoton = read np :: Int power = read pw :: Double dat <- getContents pcs <- forM (lines dat) $ \i -> do return $ (read i :: Photon) let tv = sumVectors pcs putStrLn $ show tv sumVectors :: [Photon] -> Position3 sumVectors [] = o3 sumVectors (pc:pcs) = getDirVec pc + sumVectors pcs getDirVec :: Photon -> Position3 getDirVec (_, (_, d)) = d
5回実行し、その結果のベクトルの長さを計算したところ下表の通りとなった。 もちろん理想は長さゼロだ。
| #photon | min | max | avg | err/photon | time |
|---|---|---|---|---|---|
| 10000 | 80.5 | 129.7 | 101.3 | 1.01% | 0.62 s |
| 100000 | 47.2 | 373.4 | 216.1 | 0.22% | 5.70 s |
| 1000000 | 531.8 | 1232.8 | 792.3 | 0.079% | 57.18 s |
これが統計的にみて「十分に無作為」と言えるのかどうか正直微妙だ。 しかし今の所、何か改善できるのか、どうしたらできるのかよくわからないので とりあえず放置しよう。
一方、生成されるフォトンの波長毎の数は無作為か? 測ってみたところ下記のようになった。（フォトン数10万個で5回試行）
| light | red | green | blue |
|---|---|---|---|
| 白色(2W) | 33.3% | 33.3% | 33.3% |
| 白色(2W)+橙色(1W) | 44.4% | 33.3% | 22.2% |
橙色の光源は波長の割合を赤2:緑1:青0とした。 こちらは綺麗に無作為な生成ができていると言える。まあ、波長については 一様乱数で0-1の値を生成してそれを各色の割合で直接選択しているのだから うまくいって当たり前か。
まとめ
とりあえず、フォトンを生成することはできた。次回は各フォトンを追跡して フォトンマップを作ることにしよう。
※ これまでのソースはここ。 ただし、記事より先に進んでいる場合があるので注意。
