LoginSignup
1
1

More than 5 years have passed since last update.

@eijian(Eijian)

レイトレーシング(7): 光線追跡処理の大枠

Last updated at Posted at 2017-06-06

(転載その7)

忙しさにかまけて更新をサボってしまった。。。気を取り直して再開しよう。前回まででフォトンマップが出来上がったので、今回からはそれを使った光線追跡のプログラムに取り組む。

処理の大枠

フォトンマップを使う以外は古典的な光線追跡のアルゴリズムである。一番大枠の処理は次の通り。

  1. 仮想カメラを用意（定義）する（カメラ位置、スクリーン、解像度など）。
  2. フォトンマップを読み込む。
  3. スクリーンの各画素について光線追跡を行い輝度を得る。
  4. 輝度値を画像ファイルに書き出す。

3が本プログラムのメインとなるが、今回は外堀を埋めよう。1、2、4を定義する。3も画素をループするところまで。

仮想カメラの用意

1について必要なことを書き出してみる。流儀などがあるだろうが、自分が作るときには以下の要素を用いている。

  • 視線
    • 視点
    • 視線方向: スクリーンの中心方向
    • 垂直方向: スクリーンの「上」方向を決めるための情報
  • スクリーン
    • フォーカス: スクリーンまでの距離としている
    • 解像度(x, y)

なお、スクリーンの中心点は視点から視線方向にフォーカス分の距離を進んだ地点とし、そこから縦横±1.0の領域と固定している。この正方形を解像度の設定値で細かく分割してそれぞれに光線追跡するわけだ。なお、スクリーンは視線方向に垂直な平面である。

figure-20150801.png

これらの情報は、ちゃんとしたプログラムにするときには設定ファイルなどから読み込んでいろいろ変更できるようにするのだろうが、とにかく「フォトンマッピング」法の効果を見たいので決め打ち＝ソース中に直書きしてしまう。具体的には前回示したフォトンマップの視覚化と同じ設定を使う(あとで画像を比較できるし)。今回使う設定は以下だ。

Screen.hs
eyepos = initPos 0 2 0  -- 視点
eyedir = ez3            -- 視線方向 (0, 0, 1)=Z軸
upper  = ey3            -- 上方向を表すベクトル (0, 1, 0)=Y軸  
focus  = 1.0 :: Double  -- フォーカス(スクリーンまでの距離)
xres   = 256 :: Int     -- 解像度(横方向)
yres   = 256 :: Int     -- 解像度(縦方向)

ただ、これだけでは光線追跡処理には不十分なので、いくつかの定義を追加する。最初の画素の位置や次の画素の場所（縦横方向）などを決めてやらないといけないからあらかじめ計算しておく。

stepx = 2.0 / fromIntegral xres :: Double  -- 画素の大きさ(横)
stepy = 2.0 / fromIntegral yres :: Double  -- 画素の大きさ(縦)
eex = ex3          -- スクリーンの横方向を表す方向ベクトル(大きさ1)
eey = negate ey3   -- スクリーンの縦方向を表す方向ベクトル(大きさ1)
origin = focus *> eyedir
  + ((-1.0 + 0.5 * stepx) *> eex)
  - (( 1.0 - 0.5 * stepy) *> eey)  -- スクリーンの左上の点(最初の画素)

これらを使えば、位置(X, Y)の画素pointは

point = origin + X * stepx * eex + Y * stepy * eey

で表せるわけだ。

フォトンマップの読み込み

さてほとんど忘れかけていたが、前回作ったフォトンマップのフォーマットを思い出してみよう。次のような形だった。フォトンキャッシュ情報はPhotonCache型である。

フォトン数(改行)
フォトン1個のエネルギー(改行)
1番目のフォトンキャッシュ情報(改行)
2番目のフォトンキャッシュ情報(改行)
   : 
n番目のフォトンキャッシュ情報(改行)

よって読み込みはこの逆をしてやればよい。最終的なデータ構造にはk-d treeを使う。これは例の本(フォトンマッピング)にk-d treeが良いと書いてあったのでそのまま採用しただけだ。幸い、haskellにはk-d treeのライブラリがいくつかあるので自作する必要はない。実はここでちょっとつまづいた。ライブラリによって、処理性能（速度）が段違いなのだ。ここでは最終的に採用したkdtを元に進める。最初使ったものは遅すぎてとても実用にならなかったので（多分、処理速度が数十倍違う）。インストールはcabalから。

cabal install kdt

kdtで使うため、読み込んだデータを変換してPhotonInfo型のListにしよう。PhotonInfo型の定義、変換用関数は以下のとおり（PhotonCache型を使えばいいのだが、フォトンの入射方向を逆転させる必要があるのでこれはあとで取り組む）。また、k-d treeのデータ構造へ変換するときに必要になるので、フォトンの位置を取り出す関数も定義しておく。

Optics.hs
data PhotonInfo = PhotonInfo Wavelength Position3 Direction3 deriving (Show, Eq)

convertToInfo :: PhotonCache -> PhotonInfo
convertToInfo (wl, (rp, rd)) = PhotonInfo wl rp (negate rd)

infoToPointList :: PhotonInfo -> [Double]
infoToPointList (PhotonInfo _ (Vector3 x y z) _) = [x, y, z]

あとは全フォトン情報を読み出してk-d treeに変換すれば良い。なお、フォトンマップのデータは標準入力で与えることにする。

各画素の輝度値を求める

具体的な処理は次回以降に回すとして、今回は大枠なので適当な輝度値をでっち上げ、次の処理でファイルに書き出せるようにしよう。入力は各画素の位置、出力は輝度値とする。輝度値は、各周波数の値を持った型を定義することにしよう。Radianceとする。

Optics.hs
data Radiance = Radiance Double Double Double deriving (Read, Show)

輝度値も光線追跡の最中に頻繁に演算する。なのでVectorとほとんど同じように加減などの演算関数を定義することになる。が、まずは型だけ。各画素の位置はxとyの組みで表そう。そうすると、レイトレーシングのメイン処理は大体次のような形になるだろう。全画素の組はあらかじめ作っておき、mapで各組を輝度値に変換する。

scrmap = [(y, x) | y <- [0..(yres - 1)], x <- [0..(xres - 1)]]  
  :
image = map (traceScreen ??) srcmap

traceScreen :: <<必要な引数>> -> (Int, Int) -> Radiance
traceScreen ?? scrmap = ...
  :

ひとまず、RGB各周波数の輝度値が 0.05 のRadianceのリストを返すようにしておこう。あとでわかるが、0.05にしておくと灰色になるはずだ。

輝度値の書き出し

輝度値を画像ファイルへ書き出すため、各周波数成分を256段階の整数に変換したい。次にような関数を用意すれば良いだろう。

clip :: Double
clip = 0.1
gamma = 1.0 / 2.2
rgbmax = 255.0

radianceToRgb :: Double -> Int
radianceToRgb d = floor (r * rgbmax)
  where
    d' = d / clip
    r  = (if d' > 1.0 then 1.0 else d') ** gamma

なお、clipは輝度値と最大値"255"(=rgbmax)の対応付けのための閾値。0.1という設定は、0.1[W]以上のエネルギーを持っていたら最大値とするということ。この値は適当なので、出来上がる画像の明暗を見て適当に変更するつもりだ。gammaはガンマ補正のための設定値で、一般的(?)な $\gamma = 2.2$ を用いる。

出力する画像データだが、単純なPPM形式を用いよう。フォーマットはググればすぐわかる。これを標準出力へ吐き出すようにする。以上を踏まえると、書き出し処理の関数は次のようなものになるだろう。

outputImage :: [Radiance] -> IO ()
outputImage rs = do
  putStrLn "P3"
  putStrLn "## test"
  putStrLn (show xres ++ " " ++ show yres)
  putStrLn "255"
  forM_ rs $ \i -> do
    putStrLn convertOneCell i

各画素の輝度値のリストを引数にして書きだすだけだ。なおconvertOneCellは一画素の結果を文字列に直す関数。

まとめ

上記をまとめ、ソースにしてみた。

(RT.hs)

module Main where

import System.IO
import Control.Monad
import Data.KdTree.Static
import NumericPrelude

import Ray.Algebra
import Ray.Geometry
import Ray.Physics
import Ray.Optics

-- PARAMETERS --
-- for camera
eyepos = initPos 0 2 0
eyedir = ez3
upper = ey3
focus = 1.0 :: Double
xres = 256 :: Int
yres = 256 :: Int

stepx = 2.0 / fromIntegral xres :: Double
stepy = 2.0 / fromIntegral yres :: Double
eex = ex3
eey = negate ey3
origin = focus *> eyedir
  + ((-1.0 + 0.5 * stepx) *> eex)
  - (( 1.0 - 0.5 * stepy) *> eey)

scrmap = [(y, x) | y <- [0..(yres - 1)], x <- [0..(xres - 1)]]  

-- FUNCTIONS --

main :: IO ()
main = do
  (power, photonmap) <- readMap
  let image = map traceScreen scrmap
  outputImage image
--
readMap :: IO (Double, KdTree Double PhotonInfo)
readMap = do
  np' <- getLine
  pw' <- getLine
  let np = read np' :: Int
  let pw = read pw' :: Double
  pcs <- forM ([1..np]) $ \i -> do
    l <- getLine
    return $ (read l :: PhotonCache)
  let pmap = build infoToPointList (map convertToInfo pcs)
  return (pw, pmap)

--
traceScreen :: (Int, Int) -> Radiance
traceScreen (y, x) = Radiance 0.05 0.05 0.05
--
outputImage :: [Radiance] -> IO ()
outputImage rs = do
  putStrLn "P3"
  putStrLn "## test"
  putStrLn (show xres ++ " " ++ show yres)
  putStrLn "255"
  forM_ rs $ \i -> do
    putStrLn $ convertOneCell i

convertOneCell :: Radiance -> [Char]
convertOneCell (Radiance r g b) =
  (show $ radianceToRgb r) ++ " " ++
  (show $ radianceToRgb g) ++ " " ++
  (show $ radianceToRgb b)
Physics.hs

data Wavelength = Rad | Green | Blue deriving (Show, Read, Enum, Eq)
Optics.hs

data Radiance = Radiance Double Double Double deriving (Read, Show)

clip = 0.1 :: Double
gamma = 1.0 / 2.2
rgbmax = 255.0

radianceToRgb :: Double -> Int
radianceToRgb d = floor (r * rgbmax)
  where
    d' = d / clip
    r  = (if d' > 1.0 then 1.0 else d') ** gamma


type PhotonCache = (Wavelength, Ray)
data PhotonInfo = PhotonInfo Wavelength Position3 Direction3
  deriving (Show, Eq)

infoToPointList :: PhotonInfo -> [Double]
infoToPointList (PhotonInfo _ (Vector3 x y z) _) = [x, y, z]

convertToInfo :: PhotonCache -> PhotonInfo
convertToInfo (wl, (rp, rd)) = PhotonInfo wl rp (negate rd)

コンパイルして実行してみる。

$ ghc -o rt0 RT.hs
$ rt0 < photonmap > out1.ppm
$ convert out1.ppm out1.png   <= ImageMagickを使う

結果は以下のとおり。当たり前だが面白くもなんともない・・・。

out1.png

ということで、次回は光線追跡の処理を考えていこう。

1
1
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
1
1