レイトレーシング(6): やっとフォトン追跡
シーン作成
前回までで、物体を定義するための準備ができた。それを用いて簡単なサンプルシーンを 作ってみる。閉じた空間でないと、せっかく放射したフォトンが無限遠に行ってしまって無駄に なるので、箱型の空間を用意しよう。これで平面が6つ。それだけではつまらないので、 1つだけ球面を置こう。その定義を以下に示す。
m_ball = Material 0.8 0.3 0.3 m_wall = Material 0.8 0.8 0.8 m_ceil = Material 0.4 0.2 0.02 m_flor = Material 0.8 0.6 0.4 wall_bt = initObject (Plain ey3 0) m_flor -- bottom wall_tp = initObject (Plain (negate ey3) 4) m_ceil wall_rt = initObject (Plain (negate ex3) 2) m_wall wall_lt = initObject (Plain ex3 2) m_wall wall_bk = initObject (Plain ez3 1) m_wall wall_ft = initObject (Plain (negate ez3) 5) m_wall ball1 = initObject (Sphere (initPos 1 0.8 3) 0.8) m_ball objs = [wall_bt, wall_tp, wall_rt, wall_lt, wall_bk, wall_ft, ball1]
横幅と高さが4[m]、奥行き6[m]の部屋に直径1.6[m]の球が置いてあるという想定だ。 壁は白、床はオレンジ、球は少し薄めの赤とした。まだ描画していないので色は よくわからないけど。
最終的にobjs
が物体のリストである。これを使ってフォトンを追跡しよう。
フォトンの追跡
前回
示したtracePhotons
関数を定義しよう。といっても、物体のリストとフォトンのリストを
引数にするだけだ。
tracePhotons :: [Object] -> [Photon] -> IO [PhotonCache] tracePhotons objs photons = do pcs <- mapM (tracePhoton objs) photons return $ concat pcs
うーん、何のことはない。フォトンの数だけtracePhoton
(複数形でないことに注意)を
呼び出して、結果をフラットなリストにして返しているだけだ。tracePhoton
が具体的な
追跡処理なのだが、実行することは次のようになるだろう。
これをコードにしたら次のようになった。
tracePhoton :: [Object] -> Photon -> IO [PhotonCache] tracePhoton os (wl, r) = do let iss = filter (\x -> fst x > nearly0) (concat $ map (calcDistance r) os) (t, s) = head $ sortBy (comparing fst) iss return [(wl, initRay (target t r) (getDir r))] calcDistance :: Ray -> Object -> [(Double, Object)] calcDistance r o@(Object s m) = zip ts (replicate (length ts) o) where ts = distance r s
しつこいようだが"バージョン1"ではフォトンの反射・屈折は扱わないので、上記コードは
若干無駄である。例えばtracePhoton
の返値はこの場合1個しかないので
IO [PhotonCache]
とリストにする必要はない。しかし今後、反射・屈折に対応したら
1個のフォトンの追跡で複数のフォトンキャッシュ情報が得られるようになる(ハズ)だから
今のうちにその対応をしておこうということだ。また内側で使われているcalcDistance
についても、今必要なのは交点の位置ベクトルだけなのでtarget t r
を呼んで距離を
位置ベクトルに変換するだけでよいハズだ。しかし、反射・屈折を扱うには交点での材質の
情報が必要だ。また反射・屈折方向を求めるためには衝突した物体の「形」も要る。
ということで、あえて衝突した位置までの距離と衝突した物体をセットにして返している
わけだ。
フォトンキャッシュ情報の書き出し
さあ、一連の処理の最後だ。各フォトンを追跡して得られたフォトンキャッシュ情報を 標準出力に書き出そう。これをファイルに書き込めばフォトンマップのデータファイルとなる。 ファイル形式は、
フォトン数(改行) フォトン1個のエネルギー(改行) 1番目のフォトンキャッシュ情報(改行) 2番目のフォトンキャッシュ情報(改行) : n番目のフォトンキャッシュ情報(改行)
とする。以前示したメインルーチンより少しシェイプアップして次のようになった。
putStrLn $ show nphoton putStrLn $ show power forM_ (photoncaches) $ \i -> do putStrLn $ show i
やっとできた!これまでのプログラムをまとめ、コンパイルして実行だ。 結果はかなり大きなファイルになる。10万個のフォトンで約130MBぐらい。そうか、 反射・屈折なしでこれぐらいということは、それを実装すると簡単に1GBに到達する。 まあ、ファイル形式の再考はあとにしよう。
実行結果
ここまでくると、フォトンマップがどんな結果になったのかどうしても見てみたい。 もちろん、ちゃんと動いているのか気になる所であるし。フォトンキャッシュ情報にある 交点座標を、適当な二次元平面に投影させてフォトンの色も踏まえて描き出してみた結果が 以下だ。
おお、なんとなくそれっぽいかも。ちゃんと部屋になって、右下に球も見える。 黒いところは球の影でフォトンが届かないところだ。
物は試しに、複数の光源に変えて実行してみる。天井の中央ではなく、左右奥に白と橙色の 2つの光源だ。ただ、全体の光強度は同じ、フォトン数も同じだ。
なかなかよい。とりあえず、フォトンマップは完成した、ということにしよう。
メモリ使用量の改善
・・・いや、終わらなかった・・・。
フォトン数を100万個にして処理時間を計測中、たまたま眺めていたtopコマンドの 表示に驚いた。処理中の最大メモリ使用量が約1.5GBだったのだ! 改めて作成したメインルーチンを見てみよう。
main = do (power, photons) <- generatePhotons nphoton lgts photoncaches <- tracePhotons objs photons putStrLn $ show nphoton putStrLn $ show power forM_ (photoncaches) $ \i -> do putStrLn $ show i 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)
フォトン数に応じていくつも巨大なリストを作っている。フォトンを生成するのに まず、n個の光源のリストを生成し、そのリスト全体に処理をかけている。 さらにそのフォトンのリストを処理してフォトンキャッシュのリストを生成、 最後にリスト全体を表示している。とにかく大きなリストをいくつも作っている。 全フォトンを並列に処理しているのでフォトン数が増えるほどメモリを食うのも仕方なしか。
ここで、フォトンを1個だけ追跡する処理を改めて考えてみよう。
最後にフォトンキャッシュを出力してしまえば、何も情報を保存しておく必要が ないので、並列ではなく直列に1個ずつ処理を終わらせていけばいいということか。
ような「二重ループ」にすればよいだろう。コードにしてみたのが以下。
main :: IO () main = do putStrLn $ show nphoton let power = (sum $ map flux lgts) / (fromIntegral nphoton) ns = map (calcN power) lgts putStrLn $ show power zipWithM_ outputPhotonCaches ns lgts calcN :: Double -> Light -> Int calcN power light = round (flux light / power) outputPhotonCaches :: Int -> Light -> IO () outputPhotonCaches n lgt = mapM_ (outputPhotonCache lgt) [1..n] outputPhotonCache :: Light -> Int -> IO () outputPhotonCache lgt _ = generatePhoton lgt >>= tracePhoton objs >>= mapM_ (putStrLn.show)
main
の中で求めているns
が、各光源から放出されるフォトン数n(i)のリスト。
最終行のzipWithM_
でそれぞれの光源について処理している(外側のループ)。
outputPhotonCahces
ではn(i)回フォトンキャッシュを計算して出力するように
している(内側のループ)。ちなみに2つ目の引数は使わないので捨てている。
このように「同じ処理をn回実行する」というのはどう書くのがスマート
なんだろうね。
`outputPhotnCache'(単数形)が実際に一個のフォトンを生成から追跡、最後の 書き出しまでやっている関数だ。
これをコンパイルして実行してみると、出力結果は同じだが実行時のメモリ使用量が 「数十kB」まで激減した!アルゴリズムは大事だなぁ。(というか最初考えた アルゴリズムがあまりにもイケていなかっただけか・・・)
ようやく完成した。