在随机光线追踪器的上下文中,我想将 MC 集成(路径追踪、双向路径追踪)与样本生成(均匀随机、分层、泊松、大都会......)分离。其中大部分已经实现,但是使用起来很乏味。所以我放弃了它并尝试通过将采样计算分成两个阶段来构建更好的东西:SampleGen
你可以使用mk1d
and函数请求一个随机值,然后由采样算法mk2d
提供实际的 s 。Float
可以检查这些值以进行SampleRun
实际计算。这是一些带有分层采样器有趣部分的代码,它的使用:
{-# LANGUAGE GeneralizedNewtypeDeriving #-}
import Control.Applicative
import Control.Monad.State.Strict
import Control.Monad.Primitive
import System.Random.MWC as MWC
-- allows to construct sampled computations
newtype SampleGen s m a = SampleGen (StateT s m a)
deriving ( Functor, Applicative, Monad
, MonadState s, MonadTrans )
-- allows to evaluate sampled computations constructed in SampleGen
newtype SampleRun s m a = SampleRun (StateT s m a)
deriving ( Functor, Applicative, Monad
, MonadState s )
-- a sampled computation, parametrized over the generator's state g,
-- the evaluator's state r, the underlying monad m and the result
-- type a
type Sampled g r m a = SampleGen g m (SampleRun r m a)
----------------------
-- Stratified Sampling
----------------------
-- | we just count the number of requested 1D samples
type StratGen = Int
-- | the pre-computed values and a RNG for additional ones
type StratRun m = ([Float], Gen (PrimState m))
-- | specialization of Sampled for stratified sampling
type Stratified m a = Sampled StratGen (StratRun m) m a
-- | gives a sampled value in [0..1), this is kind
-- of the "prime" value, upon which all computations
-- are built
mk1d :: PrimMonad m => Stratified m Float
mk1d = do
n1d <- get
put $ n1d + 1
return $ SampleRun $ do
fs <- gets fst
if length fs > n1d
then return (fs !! n1d)
else gets snd >>= lift . MWC.uniform
-- | gives a pair of stratified values, should really also
-- be a "prime" value, but here we just construct them
-- from two 1D samples for fun
mk2d :: (Functor m, PrimMonad m) => Stratified m (Float, Float)
mk2d = mk1d >>= \f1 -> mk1d >>= \f2 ->
return $ (,) <$> f1 <*> f2
-- | evaluates a stratified computation
runStratified
:: (PrimMonad m)
=> Int -- ^ number of samples
-> Stratified m a -- ^ computation to evaluate
-> m [a] -- ^ the values produced, a list of nsamples values
runStratified nsamples (SampleGen c) = do
(SampleRun x, n1d) <- runStateT c 0
-- let's just pretend I'd use n1d to actually
-- compute stratified samples
gen <- MWC.create
replicateM nsamples $ evalStateT x ([{- samples would go here #-}], gen)
-- estimate Pi by Monte Carlo sampling
-- mcPi :: (Functor m, PrimMonad m) => Sampled g r m Float
mcPi :: (Functor m, PrimMonad m) => Stratified m Float
mcPi = do
v <- mk2d
return $ v >>= \(x, y) -> return $ if x * x + y * y < 1 then 4 else 0
main :: IO ()
main = do
vs <- runStratified 10000 mcPi :: IO [Float]
print $ sum vs / fromIntegral (length vs)
这里缺少的部分是,在它的当前形式中,mcPi
函数具有类型
mcPi :: (Functor m, PrimMonad m) => Stratified m Float
虽然它真的应该是这样的
mcPi :: (Functor m, PrimMonad m) => Sampled g r m Float
承认,上面的四个类型参数Sampled
并不是很漂亮,但至少这样的东西会很有用。总之,我正在寻找一些允许表达计算的东西,比如mcPi
独立于采样算法,例如:
- 一个均匀随机采样器在相位中不需要保持任何状态,在相位
SampleGen
中只需要一个RNGSampleRun
- 分层和泊松圆盘采样器(可能还有其他)都跟踪所需的一维和二维样本的数量并将它们预先计算为一个向量,并且它们将被允许共享一个
SampleGen
和SampleRun
实现,只是在两者之间发生的事情有所不同SampleGen
和SampleRun
(如何实际填充向量) - Metropolis 采样器将在其阶段使用惰性样本生成技术
SampleRun
我想使用 GHC 来编译它,所以扩展喜欢MultiParamTypeClasses
并且TypeFamilies
对我来说没问题,但我没有想出任何远程可用的东西。