Haskellでいってみよう

日曜プログラマにも満たないレベルでもHaskellで何かソフトウェアを作りたい!

レイトレーシング(11): 画質向上策2点

今回は寄り道

前回から半年も経ってしまった orz

前回で完全拡散面の反射まで実装した。当然次は鏡面反射、と思われたら 申し訳ないところだが、反射(と屈折)は手間がかかるので後回し。 今回は寄り道として生成画像の画質改善に取り組む。

1. 直接光は古典レイトレーシング

フォトンマッピングでは、放射するフォトン数に限界があり、多数のフォトン から輝度推定しても"まだら模様"のようになってしまう。また、点光源では本来 影の縁はくっきりするものだが、残念ながらぼやけてしまう。

例の本(フォトンマッピング、Herik Wann Jensen、オーム社)では、 実践的なアルゴリズムとして「直接光はフォトンマップから推定するのではなく 古典レイトレーシングで」と書いてある(ようだ)。

確かにランダムなフォトンから推定するより数式できっちり計算できる 古典レイトレーシングを用いたほうが画質的には有利だろう。

前回までで、フォトンマッピング法の画像と比較するため古典レイトレーシングも 実装してきた。これをそのまま流用すればよい。基本的な考え方は次の通り。

実装

ではまずフォトン追跡の実装から。

tracePhoton :: [Object] -> Int -> Photon -> IO [PhotonCache]
tracePhoton os l (wl, r) = do
  let is = calcIntersection r os
  if is == Nothing
    then return []
    else do
      let (p, n, m) = fromJust is
      i <- russianRoulette wl [reflectance m]
      pcs <- if i > 0
        then reflect p n os l wl
        else return []
      if (useClassicForDirect == False || l > 0) && diffuseness m > 0.0
        then return $ ((wl, initRay p (getDir r)) : pcs)
        else return pcs

tracePhotonに$l$という第二引数を設けている。これは追跡の深さを 表すもので、物体に当たって反射するたびに1ずつ増やして反射方向の追跡を するようにした。つまり$l=0$ならまだ反射していないから「光源からの 直接光」である。

  if (useClassicForDirect == False || l > 0) && diffuseness m > 0.0

これがその判定をしているところ。

  • useClassicForDirectは、光源からの直接光を古典レイトレーシングで計算 するかどうかのフラグ。Falseはこれまでどおりフォトンマップから推定する という意味。
  • l > 0は前述の通り一回以上反射したフォトンを意味する。
  • diffuseness m > 0.0フォトンが衝突した物体表面が、少しでも 拡散反射の要素を持っている=フォトンを記録する必要があるという意味。

これらのいずれにも合致しない場合は、そのフォトンの衝突点ではフォトンマップに 記録しない。

次にレイトレーシング部分を見てみよう。

traceRay :: Int -> Double -> KT.KdTree Double PhotonInfo -> [Object]
         -> [Light] -> Ray -> IO Radiance
traceRay 10 _ _ _ _ _ = return radiance0
traceRay l pw pmap objs lgts r
  | is == Nothing = return radiance0
  | otherwise     = return (em + di + ii)
  where
    is = calcIntersection r objs
    (p, n, m) = fromJust is
    em = sr_half *> emittance m
    di = if useClassicForDirect
      then brdf m $ foldl (+) radiance0 $ map (getRadianceFromLight objs p n) lg ts
      else radiance0
    ii = estimateRadiance pw pmap (p, n, m)

where節の後半でuseClassicForDirectというフラグをチェックしている。 ここの仕組みはこうだ。

  • useClassicForDirectが真の場合は、物体上の点に注ぐ直接光を古典 レイトレーシングで計算する。getRadianceFromLight関数がそれ。光源の 数だけ繰り返している。偽の場合は、古典では輝度を計算しないので値ゼロ (=radiance0)としている。
  • その点には間接光も入ってくるので、その輝度は従来通りフォトンマップから 推定する。(ifの下のestimateRadiance)

なおフォトンマップから求められる輝度は、古典レイトレーシングを使う場合には 「直接光のフォトンを記録しない」ようにしてあるので、ちゃんと間接的に注ぐ 光だけで輝度推定される。古典を使わない場合は直接光のフォトンも記録されて いるので従来と同じ結果となる。

画像比較

では生成された画像を比べてみる。まずは点光源の場合。

f:id:eijian:20180108095118p:plain

効果は一目瞭然だろう。直接光を古典で計算した場合(図2)は奥の壁、床、球の 上面が滲みなく滑らかに描画されている。直接光の当たらない天井や影の部分は 間接光が支配的なので致し方ないが。注目は球の影。従来はどうしても縁がぼやけて しまったがそこがくっきりしている。点光源はこうでないと。次に平面光源の場合。

f:id:eijian:20180108095144p:plain

こちらも同様にスムーズだ。最後に平行光源の場合。

f:id:eijian:20180108095153p:plain

図5は壁面の直接光の当たる境界がぼやけてしまっているが、古典を使うと くっきりしたエッジになっている。本来こうあるべきだ。 (こちらの画像の場合天窓からの光以外は全て間接光なので、ほとんどマダラが 残ったままだが、ここは我慢)

処理時間も比べてみよう。直接光を古典レイトレーシングに任せることで 直接光のフォトンを記録しないと書いた。つまりフォトンマップが小さく なるということだ。結果は次の表のとおり。

画像 直接光 フォトン 処理時間(MM:SS)
1 点光源 フォトンマップ 318,164 6:50
2 点光源 古典レイトレ 117,852 1:22
3 面光源 フォトンマップ 288,976 6:14
4 面光源 古典レイトレ 89,420 1:09
5 平行光源 フォトンマップ 321,142 23:07
6 平行光源 古典レイトレ 120,841 2:01

どうやらフォトンマップが小さいと輝度推定時にフォトン収集の時間が短くて 済むようで、圧倒的に処理時間が短くなることがわかった! 特に平行光源の例では1/10に短縮できている。なぜこれほど短縮できるのか (もしくはフォトンマップだけだと時間がかかるのか)はよくわからないが。

2. えせアンチエリアシング

もう一つの画質改善ネタであるアンチエリアシングに移ろう。 アンチエリアシングの詳しい説明はWebを漁れば山ほど出てくるのでそちらを 見ていただくとして、今回実装した仕組みを説明しよう。

仕組み

レイトレーシングは画面を格子状に分割し、その分割された小さい四角 (ピクセル)の中心に向かって光を逆向きに追跡する(追跡レイとする)。 各ピクセルの輝度値がある閾値を超えたら、そのピクセルの内側のいくつかの 点を追加で追跡し、その輝度を平均する。

節タイトルで「えせ」としたのは、ちゃんとアンチエリアシングしたければ 追加する追跡レイはピクセル内の「ランダム」な場所に対して行うべきだから である。規則的に追跡レイを追加するとどうしてもエリアシングの問題が 残るらしい。しかし今回は、小難しいことを考えるのがだるかったので実装の 簡便さを優先した。。。時間があったらちゃんと実装しよう。

まず、すべてのピクセルで追跡レイを追加していたら時間ばかりかかって非効率 なので追加すべきかどうかをまず判断する。そのため、とあるピクセル(p)とその 周りの8つのピクセルについてそれぞれ差を求め、閾値を超えているか判定する。 1点でも閾値を超えていたら追加レイを飛ばす。(下図)

f:id:eijian:20180108095205p:plain

なお、このプログラムでは差を求める際に輝度を比較するのではなくRGB値、 すなわち「画面に出力する色」に変換してから比較している。なぜなら、 結局目に見えるところでの色の違いがある場合に追加レイを飛ばしたいから。 最初は輝度で試したが、閾値の設定がシーンの明るさ・光の量で変える必要が あり設定が難しく、実用的でないと判断した。

次に追加レイを飛ばす場所だが、上記の通り規則的にした。ピクセル内の 中心以外の4点だ。単純にピクセルを4分割してそれぞれの中心点に向かって 追跡している。

f:id:eijian:20180108095212p:plain

これら4点と最初に計算した中心点の全5点から得られた輝度を平均するのだが、 その方法には単純平均するか重み付けするかがある。結果に微妙な差が出るが 大きな問題はないと判断し、今回は単純平均する。ランダムに追加レイを飛ばす 場合は中心からの距離に差が出ると思うのでちゃんと重み付けした方が よいだろう。

実装

アンチエリアシングを実施するかどうかはメインループの最後のところ。

  :
  forM_ [0..(V.length pixels - 1)] $ \i -> do
    rgb <- smooth antiAliasing tracer pixels i
    putStrLn $ rgbToString rgb

smooth関数で処理している。引数のantiAliasingアンチエリアシングを 実施するかどうかのフラグ。smoothは次の通り。

smooth :: Bool -> (Ray -> IO Radiance) -> V.Vector Rgb -> Int -> IO Rgb
smooth False _ ims i = return (ims V.! i)
smooth True tracer ims i
  | isDifferent i ims imgOffset == False = return (ims V.! i)
  | otherwise = do
    l <- retrace tracer i
    return $ avg ((ims V.! i):l)

まずフラグが偽なら何もせず最初に求めた色を返す。真の場合、まずは追跡レイを 追加するかどうか判断する(isDifferent)。追加が必要なければやはりそのまま 色を返す。追加が必要となったら再トレースだ(retrace)。そこで得られた色を avgで平均して返せばOK。

追加が必要かどうかの判断は次の通り。

isDifferent :: Int -> V.Vector Rgb -> [Int] ->  Bool
isDifferent _ _ [] = False
isDifferent p rs (i:is)
  | p' <  0         = isDifferent p rs is
  | p' >= length rs = isDifferent p rs is
  | df              = True
  | otherwise       = isDifferent p rs is
  where
    p' = p + i
    df = diffRgb (rs V.! p) (rs V.! p')

少し説明が必要かもしれない。2番目の引数は、すでに 計算済みの画面全部の色値(RGB)の配列(Vector)である。着目しているピクセルは 第一引数pに配列内の位置として与えられる。p'が比較対象となる 周りのピクセル(の位置)は、第3引数で与えられている。 これは上の方で次のように定義されている。

imgOffset :: [Int]
imgOffset = [-xres-1, -xres, -xres+1, -1, 1, xres-1, xres, xres+1]

画面の横方向画素数xresを使い、着目点pに対するオフセットを 示しているわけだ。このようにして、周りの8個の点と比べ、閾値を 超えたら直ちにTrueを返すようにしてある。周りの点が画面をはみ出す 場合はもちろん無視だ(条件 p' < 0p' >= length rs)。

次に再トレース部分を説明しよう。

retrace :: (Ray -> IO Radiance) -> Int -> IO [Rgb]
retrace tracer p = do
  let p' = (fromIntegral (p `div` xres), fromIntegral (p `mod` xres))
  ls <- mapM tracer $ map generateRay' (map (badd p') blur)
  return $ map radianceToRgb ls

といっても大したことはない。現在の着目点pに対し、1/4ずつずれた4点を 求め、それぞれ普通に追跡して結果をリストにして返しているだけだ。 ずれた点を計算しやすくするため、blurをあらかじめ定義している。

blur :: [(Double, Double)]
blur = [(-0.25, -0.25), (-0.25,  0.25), ( 0.25, -0.25), ( 0.25,  0.25)]

生成画像の比較

では実際にアンチアリアスを効かせた画像を生成してみよう。従来の画像との 比較はつぎのとおり。

f:id:eijian:20180108095220p:plain

f:id:eijian:20180108095229p:plain

ぱっと見はよくわからないかもしれないが、特に球の上部や平面光源の縁が気持ち なだらかになったのがわかるだろうか?これがアンチエリアシングの効果だ。 ただ、いくつかのピクセルについて追加レイを飛ばすので、当然計算時間が 余計にかかる。

アンチエリアシング フォトン 処理時間(MM:SS)
89420 1:09
89420 1:52

画面の複雑さや色合いにもよるので一般にどのぐらい時間が増えるかは なんとも言えないが、この例ではほぼ2倍の計算時間がかかったことになる。 この辺は画質とのトレードオフだが、世間のCGソフトは普通アンチエリアシングが 効いているからなぁ。

まとめ

今回は多少なりとも画質を向上する策として、以下の2点の改善を実施した。

両方とも、簡便な方法を採った割にはそれなりにうまくいったと思う。 また、古典レイトレーシングを使うことで大幅な計算時間の短縮ができた ことは収穫だった。

次回はとうとう(?)鏡面反射、屈折に対応しよう。うまく行けばいわゆる 「集光模様」が得られるはずだ。


ここまでのソースはこちら